# loading packages
# devtools::install_github("thomasp85/patchwork")
pacman::p_load(tidyverse, # tidy family and related pacakges below
kableExtra,
gridExtra, # may not use this
purrr,
magrittr, # extending piping
pander, # nice tables
metafor, # package for meta-analysis
MCMCglmm, # Bayeisan mixed model package
ggbeeswarm, # making bee-swarm plots possible
plotly, # interactive plots using ggplot2
MuMIn, # multi-model inference
lme4, # lmm & glmm (models)
broom.mixed, # getting estimates from lmer + glmer objects
performance, # getting R2 from lmer + glmer objects
png, # reading png files
grid, # graphic layout manipulation
patchwork, # putting ggplots together - you need to install via devtool
here # making reading files easy
#lmerTest, # more functions for lme4
#mi, # missing data analysis
#betareg # dependence of the above
)
We have 5 custom functions named : p_to_Zr(),I2(), R2(), get_est(), get_pred(), and cont_gen(), all of which are used later (see below for their functionality) and the code is included here.
# coustm functions
#' Title: getting Zr and its sampling variance from p value
#'
#' @param data: data frame
#' @param pval: p value
#' @param N: sample size (N: the number of species ) and the degrees of freedom df = N - 2
#'
#' @return
#' @export
#'
#' @examples
p_to_Zr <- function(data, pval, N) {
# turning them into strings
pval <- data[[deparse(substitute(pval))]]
N <- data[[deparse(substitute(N))]]
# getting t values
tval <- -qt(pval, N - 2)
rval <- tval/sqrt((tval^2) + (N - 2))
# define Zr function Zr <- 0.5*(log(1 + rval) - log(1 - rval)); the same as below
# r <-tanh(Zr) # turning Zr to r
Zr <- atanh(rval)
# getting Var(Zr)
VZr <- 1/(N - 3)
# putting all together
Zrs <- tibble(rval, Zr, VZr)
data <- bind_cols(data, Zrs)
}
# coverting back Zr to r Just use 'psych' pacakge - fisherz2r(z) -
# <http://personality-project.org/r/psych/help/fisherz.html> or this will do : r
# to Zr is tanh(r)!!
# Functions for processing
# General modeling functions Functions for I2
#' Title Function to obtain total and separate I2 from multilevel-meta-analytic model
#'
#' @param model
#' @param method
#'
#' @return
#' @export
#'
#' @examples
I2 <- function(model, method = c("Wolfgang", "Shinichi")) {
## evaluate choices
method <- match.arg(method)
# Wolfgang's method
if (method == "Wolfgang") {
W <- solve(model$V)
X <- model.matrix(model)
P <- W - W %*% X %*% solve(t(X) %*% W %*% X) %*% t(X) %*% W
I2_total <- sum(model$sigma2)/(sum(model$sigma2) + (model$k - model$p)/sum(diag(P)))
I2_each <- model$sigma2/(sum(model$sigma2) + (model$k - model$p)/sum(diag(P)))
names(I2_each) = paste0("I2_", model$s.names)
# putting all together
I2s <- c(I2_total = I2_total, I2_each)
# or my way
} else {
# sigma2_v = typical sampling error variance
sigma2_v <- sum(1/model$vi) * (model$k - 1)/(sum(1/model$vi)^2 - sum((1/model$vi)^2))
I2_total <- sum(model$sigma2)/(sum(model$sigma2) + sigma2_v) #s^2_t = total variance
I2_each <- model$sigma2/(sum(model$sigma2) + sigma2_v)
names(I2_each) = paste0("I2_", model$s.names)
# putting all together
I2s <- c(I2_total = I2_total, I2_each)
}
return(I2s)
}
# test <- dataset$fit4.1[[3]] I2(test, method = 'Wolfgang') I2(test, method =
# 'Shinichi')
#' Title: R2 based on Nakagawa & Schielzeth 2013
#'
#' @param model
#'
#' @return
#' @export
#'
#' @examples
R2 <- function(model) {
warning("Conditional R2 is not meaningful and the same as marginal R2\n")
# fixed effect variance
fix <- var(as.numeric(as.vector(model$b) %*% t(as.matrix(model$X))))
# marginal
R2m <- fix/(fix + sum(model$sigma2))
R2
# Rm <- round(100*R2m, 3)
# conditional
R2c <- (fix + sum(model$sigma2) - model$sigma2[length(model$sigma2)])/(fix +
sum(model$sigma2))
R2s <- c(R2_marginal = R2m, R2_coditional = R2c)
return(R2s)
}
#' Title: the function to get estimates from rma objects (metafor)
#'
#' @param model: rma.mv object
#' @param mod: the name of a moderator
get_est <- function(model, mod = " ") {
name <- as.factor(str_replace(row.names(model$beta), mod, ""))
estimate <- as.numeric(model$beta)
lowerCL <- model$ci.lb
upperCL <- model$ci.ub
table <- tibble(name = name, estimate = estimate, lowerCL = lowerCL, upperCL = upperCL)
}
#' Title: the function to get prediction intervals (crediblity intervals) from rma objects (metafor)
#'
#' @param model: rma.mv object
#' @param mod: the name of a moderator
get_pred <- function(model, mod = " ") {
name <- as.factor(str_replace(row.names(model$beta), mod, ""))
len <- length(name)
if (len != 1) {
newdata <- matrix(NA, ncol = len, nrow = len)
for (i in 1:len) {
# getting the position of unique case from X (design matrix)
pos <- which(model$X[, i] == 1)[[1]]
newdata[, i] <- model$X[pos, ]
}
pred <- predict.rma(model, newmods = newdata)
} else {
pred <- predict.rma(model)
}
lowerPR <- pred$cr.lb
upperPR <- pred$cr.ub
table <- tibble(name = name, lowerPR = lowerPR, upperPR = upperPR)
}
# Here are links for how to do confidence regions for rma.mv regression lines
# https://www.rdocumentation.org/packages/metafor/versions/1.9-9/topics/predict.rma
# https://stackoverflow.com/questions/50804464/out-of-sample-prediction-for-rma-object-in-metafor
#' Title: Contrast name geneator
#'
#' @param name: a vector of character strings
cont_gen <- function(name) {
combination <- combn(name, 2)
name_dat <- t(combination)
names <- paste(name_dat[, 1], name_dat[, 2], sep = "-")
return(names)
}
Citations for papers describing the main methods of cophylogeny analysis, based on a Google Scholar search conducted on 4th July 2019. Although the paper describing Brooks parsimony analysis has more citations than that describing the Parafit method, relatively few citations correspond to actual cophylogenetic analyses employing the approach, versus methodological discussion.
# getting the data and formating some variables (turning chraracter vectors to
# factors)
read_csv(here("data/lit_search.csv"), na = "NA") %>% mutate_if(is.character, as.factor) %>%
kable("html") %>% kable_styling("striped", position = "left")
| Method | Paper | No. citations |
|---|---|---|
| TreeMap | Page RDM. 1994. Parallel phylogenies: reconstructing the history of host–parasite assemblages. Cladistics 10: 155–173. | 385 |
| Brooks Parsimony Analysis (BPA) | Brooks DR. 1981. Hennig’s parasitological method: a proposed solution. Systematic Zoology 30: 229–249. | 362 |
| ParaFit | Legendre P, Desdevises Y, Bazin E. 2002. A statistical test for host–parasite coevolution. Systematic Biology 51: 217–234. | 344 |
| JANE | Conow C, Fielder D, Ovadia Y, Libeskind-Hadas R. 2010. Jane: a new tool for the cophylogeny reconstruction problem. Algorithms for Molecular Biology 5: 16. | 249 |
| TREEFITTER | Ronquist F. 1995. Reconstructing the history of host–parasite associations using generalised parsimony. Cladistics 11: 73–89. | 124 |
| PACO | Balbuena JA, Míguez-Lozano R, Blasco-Costa I. 2013. PACo: a novel procrustes application to cophylogenetic analysis. PloS one. 8(4):e61048. | 103 |
| TARZAN | Merkle D, Middendorf M. 2005. Reconstruction of the cophylogenetic history of related phylogenetic trees with divergence timing information. Theory in Biosciences 123: 277–299. | 72 |
| COALA | Baudet C, Donati B, Sinaimeri B, Crescenzi P, Gautier C, Matias C, Sagot MF 2014. Cophylogeny reconstruction via an approximate Bayesian computation. Systematic Biology, 64(3), 416-431. | 23 |
| COMPONENT | Page RDM. 1993. User’s manual for component, version 2.0. London, UK: The Natural History Museum. | 13 |
Below is the dataset used for our meta-analysis, followed by explanations of 24 variables extracted from the papers included (not all variables were used for our analyses; variables which were neither ‘directly’ nor ‘indirectly’ used in our analyses are indicated by *).
The meta-analytic dataset of this study.
# getting the data and formating some variables (turning chraracter vectors to
# factors)
full_data <- read_csv(here("data/2021-09-01-source-data-dat.csv"), na = "NA") %>%
mutate_if(is.character, as.factor)
# dataset to compare the same cophylogenies between the two methods
full_pair <- read_csv(here("data/2020-08-12-paried.csv"), na = "NA") %>% mutate_if(is.character,
as.factor)
# making a scrollable table
kable(full_data, "html") %>% kable_styling("striped", position = "left") %>% scroll_box(width = "100%",
height = "500px")
| authors | year | host_tax_broad | symbiont_tax_broad | symbiont_euk | symbiosis | endo_or_ecto | mode_of_transmission_broad | mode_of_transmission_fine | Visiting_symbiont? | host_tips_linked | host_tips_linked_corrected | host_genera | total_host_symbioint_links | host_range_link_ratio | host_range_taxonomic_breadth | symbiont_tips_linked | symbiont_genera | no_randomizations | p_value | method |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Althoff_et_al_2012 | 2012 | Plant | Invert | y | Mutualist | Ecto | horizontal | autonomous | visitor | 24 | 24 | 1 | 40 | 2.00 | 1.35 | 20 | 1 | 1000 | 0.00100 | Parafit |
| Althoff_et_al_2012 | 2012 | Plant | Invert | y | Mutualist | Ecto | horizontal | autonomous | visitor | 24 | 24 | 1 | 38 | 2.24 | 1.53 | 17 | 1 | 1000 | 0.00100 | Parafit |
| Arab_et_al_2019 | 2019 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 55 | 55 | 52 | 55 | 1.00 | 1.00 | 55 | 1 | 999 | 0.00100 | Parafit |
| Ballinger_et_al_2018 | 2018 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 11 | 11 | 1 | 12 | 1.00 | 1.00 | 12 | 1 | 999 | 0.07000 | Parafit |
| Banks_et_al_2006 | 2006 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 18 | 18 | 6 | 30 | 2.00 | 1.60 | 15 | 2 | 10000 | 0.00100 | Parafit |
| Bayerlova_2009 | 2009 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 21 | 21 | 14 | 31 | 1.55 | 1.45 | 20 | 1 | 9999 | 0.00040 | Parafit |
| Bellec_et_al_2014 | 2014 | Plant | Microbe | n | Parasite | Endo | horizontal | contact | resident | 22 | 22 | 3 | 133 | 2.61 | 1.65 | 51 | 1 | 999 | 0.00100 | Parafit |
| Bruyndonckxx_et_al_2009 | 2009 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 20 | 20 | 7 | 21 | 1.91 | 2.27 | 11 | 2 | 9999 | 0.00300 | Parafit |
| Caraguel_et_al__2007 | 2007 | Microbe | Microbe | y | Mutualist | Endo | vertical | vertical | resident | 6 | 6 | 1 | 6 | 1.00 | 1.00 | 6 | 1 | 9999 | 0.00100 | Parafit |
| Carneiro_et_al_2018 | 2018 | Vert | Microbe | n | Parasite | Endo | horizontal | body fluid | resident | 26 | 26 | 21 | 26 | 1.00 | 1.00 | 26 | 1 | 100000 | 0.01000 | Parafit |
| Catanach_et_al_2018 | 2018 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 54 | 28 | 5 | 44 | 1.07 | 0.93 | 41 | 30 | 999 | 0.00100 | Parafit |
| Chen_et_al_2017 | 2017 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 50 | 50 | 11 | 50 | 1.00 | 1.00 | 50 | 1 | 999 | 0.00100 | Parafit |
| Choi_&_Thines_2015 | 2015 | Plant | Microbe | y | Parasite | Endo | horizontal | autonomous | resident | 63 | 63 | 28 | 63 | 1.00 | 1.00 | 63 | 3 | 999 | 0.00100 | Parafit |
| Conord_et_al_2008 | 2008 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 14 | 14 | 10 | 14 | 1.00 | 1.00 | 14 | 1 | 999 | 0.00900 | Parafit |
| Cornuault_et_al_2012 | 2012 | Vert | Microbe | y | Parasite | Endo | horizontal | vector | resident | 8 | 8 | 1 | 23 | 1.28 | 1.17 | 18 | 1 | 9999 | 0.03500 | Parafit |
| Cruaud_et_al_2012 | 2012 | Plant | Invert | y | Mutualist | Endo | horizontal | autonomous | resident | 200 | 200 | 1 | 200 | 1.00 | 1.00 | 200 | 20 | 9999 | 0.01000 | Parafit |
| Cui_et_al_2014 | 2014 | Vert | Microbe | n | Parasite | Endo | horizontal | body fluid | resident | 46 | 46 | 46 | NA | NA | NA | 61 | NA | 99999 | 0.23300 | Parafit |
| Deng_et_al_2013 | 2013 | Invert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 7 | 7 | 4 | 10 | 1.00 | 1.00 | 10 | 2 | 999 | 0.01602 | Parafit |
| Desdevises_et_al_2002 | 2002 | Vert | Invert | y | Parasite | Ecto | horizontal | autonomous | resident | 14 | 14 | 11 | 39 | 1.95 | 1.65 | 20 | 2 | 999 | 0.26000 | Parafit |
| Dhami_et_al_2013 | 2013 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 42 | 42 | NA | 42 | 1.00 | 1.00 | 42 | 4 | 999 | 0.00100 | Parafit |
| Dona_2018 | 2018 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 14 | 14 | 13 | 15 | 1.00 | 1.00 | 15 | 1 | 100000 | 0.01000 | Parafit |
| Dona_2018 | 2018 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 42 | 42 | 29 | 44 | 1.00 | 1.00 | 44 | 1 | 100000 | 0.01000 | Parafit |
| Dowie_et_al_2016 | 2016 | Microbe | Plant | y | Parasite | Endo | horizontal | environmental | resident | 4 | 4 | 1 | 9 | 2.25 | 1.50 | 4 | 1 | 9999 | 0.00100 | Parafit |
| Du_Toit_et_al_2013 | 2013 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 4 | 4 | 1 | 14 | 1.17 | 1.17 | 12 | 1 | 10000 | 0.88000 | Parafit |
| Endara_et_al_2017 | 2017 | Plant | Invert | y | Parasite | Ecto | horizontal | autonomous | visitor | 18 | 18 | 1 | 29 | 1.00 | 2.17 | 29 | NA | 100 | 0.70000 | Parafit |
| Endara_et_al_2017 | 2017 | Plant | Invert | y | Parasite | Ecto | horizontal | autonomous | visitor | 10 | 10 | 1 | 18 | 1.50 | 1.25 | 12 | NA | 100 | 0.90000 | Parafit |
| Endara_et_al_2017 | 2017 | Plant | Invert | y | Parasite | Ecto | horizontal | autonomous | visitor | 14 | 14 | 1 | 18 | 1.20 | 1.13 | 15 | NA | 100 | 0.74000 | Parafit |
| Endara_et_al_2018 | 2018 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 44 | 30 | 1 | 45 | 1.18 | 1.03 | 38 | NA | 9999 | 0.01500 | Parafit |
| FerrerParis_et_al_2013 | 2013 | Plant | Invert | y | Parasite | Ecto | horizontal | autonomous | visitor | 64 | 64 | NA | 112 | 2.73 | 5.20 | 41 | NA | 999 | 0.15700 | Parafit |
| FraijaFernandez_et_al_2016 | 2016 | Vert | Invert | y | Parasite | Endo | horizontal | trophic | resident | 31 | 31 | 24 | 50 | 5.56 | 3.67 | 9 | 6 | 999 | 0.00100 | Parafit |
| Garamszegi_2009 | 2009 | Vert | Microbe | y | Parasite | Endo | horizontal | vector | resident | 23 | 23 | 23 | 43 | 2.39 | 2.56 | 18 | 1 | 1000 | 0.00100 | Parafit |
| Garcia_&_Hayman_2016 | 2016 | Vert | Microbe | y | Parasite | Endo | horizontal | trophic | resident | 22 | 22 | 22 | 36 | 1.33 | 1.96 | 27 | 1 | 999 | 0.01000 | Parafit |
| Gavotte_et_al_2007 | 2007 | Microbe | Microbe | n | Parasite | Endo | both | NA | resident | 33 | 33 | 1 | 51 | 0.93 | 1.00 | 55 | 1 | 10000 | 0.13190 | Parafit |
| Goker_et_al_2011 | 2011 | Microbe | Microbe | n | Parasite | Endo | NA | NA | resident | 8 | 8 | 8 | 8 | 1.00 | 1.00 | 8 | 1 | 9999 | 0.09780 | Parafit |
| Gomard_et_al_2016 | 2016 | Vert | Microbe | n | Parasite | Endo | horizontal | NA | resident | 12 | 12 | 11 | 26 | 1.00 | 1.00 | 26 | 1 | 999 | 0.09000 | Parafit |
| Gottschling_et_al_2011 | 2011 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 43 | 43 | 40 | 78 | 1.00 | 1.00 | 78 | 30 | 9999 | 0.00010 | Parafit |
| Graca_et_al_2018 | 2018 | Vert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 9 | 9 | 7 | 18 | 1.29 | 1.29 | 14 | 1 | 10000 | 0.00010 | Parafit |
| Hall_et_al__2016 | 2016 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 37 | 37 | 18 | 37 | 1.00 | 1.00 | 37 | 1 | 10000 | 0.00100 | Parafit |
| Hall_et_al__2016 | 2016 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 20 | 20 | 9 | 20 | 1.00 | 1.00 | 20 | 1 | 10000 | 0.38700 | Parafit |
| Hammer_et_al_2010 | 2010 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 23 | 23 | 6 | 23 | 1.10 | 1.00 | 21 | 1 | 999 | 0.00100 | Parafit |
| Hammerlinck_et_al_2016 | 2016 | Invert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 11 | 11 | 2 | 11 | 1.38 | 1.13 | 8 | 2 | 999 | 0.12900 | Parafit |
| Hammerlinck_et_al_2016 | 2016 | Invert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 6 | 6 | 2 | 6 | 1.50 | 1.25 | 4 | 2 | 999 | 0.19500 | Parafit |
| Hammerlinck_et_al_2016 | 2016 | Invert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 7 | 7 | 2 | 7 | 1.00 | 1.00 | 7 | 2 | 999 | 0.04200 | Parafit |
| Hembry_et_al_2013 | 2013 | Plant | Invert | y | Mutualist | Ecto | horizontal | autonomous | visitor | 37 | 37 | 1 | 35 | 1.00 | 1.00 | 35 | 1 | 999 | 0.01700 | Parafit |
| Herrera__et_al_2016 | 2016 | Microbe | Microbe | y | Parasite | Endo | horizontal | environmental | resident | 13 | 13 | 9 | 13 | 1.00 | 1.00 | 13 | 3 | 999 | 0.00500 | Parafit |
| Hewitt_et_al_2019 | 2019 | Vert | Invert | y | Parasite | Ecto | horizontal | autonomous | resident | 178 | 178 | 79 | 495 | 7.17 | 3.30 | 69 | 35 | 999 | 0.00100 | Parafit |
| Hoglund_et_al_2003 | 2003 | Vert | Invert | y | Parasite | Endo | horizontal | trophic | resident | 8 | 8 | 8 | 10 | 2.00 | 2.40 | 5 | 1 | 10000 | 0.09800 | Parafit |
| Holzer_et_al_2018 | 2018 | Invert | Invert | y | Parasite | Endo | horizontal | environmental | resident | 23 | 23 | 22 | 39 | 1.00 | 1.00 | 39 | 21 | 1000 | 0.00100 | Parafit |
| Holzer_et_al_2018 | 2018 | Vert | Invert | y | Parasite | Endo | horizontal | environmental | resident | 69 | 69 | 62 | 101 | 1.00 | 1.00 | 101 | 15 | 1000 | 0.00100 | Parafit |
| Holzer_et_al_2018 | 2018 | Vert | Invert | y | Parasite | Endo | horizontal | environmental | resident | 69 | 69 | 58 | 75 | 1.00 | 1.00 | 75 | 21 | 1000 | 0.00100 | Parafit |
| Hughes_et_al_2007 | 2007 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 18 | 18 | 6 | 18 | 1.06 | 1.00 | 17 | 1 | 999 | 0.00010 | Parafit |
| Huyse_&_Volckaert_2005 | 2005 | Vert | Invert | y | Parasite | Ecto | horizontal | autonomous | resident | 8 | 8 | 3 | 22 | 1.29 | 1.29 | 17 | 1 | 999 | 0.09500 | Parafit |
| Irwin_et_al_2012 | 2012 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 21 | 21 | 14 | 31 | 1.55 | 1.85 | 20 | NA | 9999 | 0.00040 | Parafit |
| Jenkins_et_al_2012 | 2012 | Vert | Microbe | y | Parasite | Endo | horizontal | autonomous | resident | 52 | 52 | 40 | 138 | 1.45 | NA | 95 | 1 | 999 | 0.00100 | Parafit |
| Jesovnik_et_al_2017 | 2017 | Invert | Microbe | y | Mutualist | Ecto | both | environmental | resident | 32 | 6 | 1 | 6 | 1.00 | 1.00 | 6 | 1 | 999 | 0.02900 | Parafit |
| Jousselin_et_al_2008 | 2008 | Plant | Invert | y | Mutualist | Endo | horizontal | autonomous | resident | 15 | 15 | 1 | 15 | 1.07 | 1.07 | 14 | 6 | 9999 | 0.00300 | Parafit |
| Jousselin_et_al_2008 | 2008 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 13 | 13 | 1 | 13 | 1.00 | 1.00 | 13 | 1 | 9999 | 0.02000 | Parafit |
| Jousselin_et_al_2008 | 2008 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 16 | 16 | 1 | 18 | 1.06 | 1.06 | 17 | 2 | 9999 | 0.00100 | Parafit |
| Jousselin_et_al_2008 | 2008 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 13 | 13 | 1 | 14 | 1.00 | 1.00 | 14 | 1 | 9999 | 0.00300 | Parafit |
| Jousselin_et_al_2009 | 2009 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 55 | 22 | 1 | 22 | 1.00 | 1.00 | 22 | 1 | 9999 | 0.00100 | Parafit |
| Kaltenpoth_et_al_2014 | 2014 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 39 | 39 | 2 | 41 | 1.00 | 1.00 | 41 | 1 | 1000 | 0.00100 | Parafit |
| Kawakita_&_Kato_2009 | 2009 | Plant | Invert | y | Mutualist | Ecto | horizontal | autonomous | visitor | 10 | 10 | 7 | 10 | 1.00 | 1.00 | 10 | 1 | 100 | 0.39000 | Parafit |
| Kawakita_et_al_2004 | 2004 | Plant | Invert | y | Mutualist | Ecto | horizontal | autonomous | visitor | 18 | 18 | 1 | 18 | 1.00 | 1.00 | 18 | 1 | 999 | 0.00500 | Parafit |
| Kawazoe_et_al_2008 | 2008 | Invert | Invert | y | Mutualist | Ecto | vertical | vertical | resident | 4 | 4 | 1 | 5 | 1.25 | 1.25 | 4 | 1 | 9999 | 0.00010 | Parafit |
| Kolsch_&__Pedersen_2010 | 2010 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 41 | 41 | 7 | 41 | 1.00 | 1.00 | 41 | 1 | 9999 | 0.00100 | Parafit |
| Krasnov_&_Shenbrot_2013 | 2013 | Vert | Invert | y | Parasite | Ecto | horizontal | contact | resident | 21 | 21 | 8 | 62 | 3.26 | 4.37 | 19 | 7 | 999 | 0.16000 | Parafit |
| Krumbholz_et_al_2009 | 2009 | Vert | Microbe | n | Parasite | Endo | horizontal | NA | resident | 13 | 13 | 13 | 18 | 1.00 | 1.00 | 18 | 1 | 999999 | 0.49460 | Parafit |
| Ku_&_Hu_2014 | 2014 | Plant | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 11 | 11 | 1 | 11 | 1.00 | 1.00 | 11 | 1 | 9999 | 0.02000 | Parafit |
| Lanterbecq_et_al_2010 | 2010 | Invert | Invert | y | Parasite | Endo/Ecto | horizontal | contact | resident | 16 | 16 | 12 | 16 | 1.00 | 1.00 | 16 | 5 | 1000 | 0.01300 | Parafit |
| Latinne_et_al_2018 | 2018 | Vert | Microbe | y | Parasite | Endo | horizontal | enviromental | resident | 15 | 15 | 6 | 19 | 3.17 | 2.00 | 6 | 1 | 999 | 0.00900 | Parafit |
| Lauber_et_al_2017 | 2017 | Vert | Microbe | n | Parasite | Endo | horizontal | autonomous | resident | 30 | 30 | 28 | 30 | 0.88 | 1.00 | 34 | 8 | 10000 | 0.00010 | Parafit |
| Lauron_et_al_2015 | 2015 | Vert | Microbe | y | Parasite | Endo | horizontal | vector | resident | 18 | 18 | 8 | 83 | 1.80 | NA | 46 | 1 | 1000 | 0.66000 | Parafit |
| Lei_&_Olival_2014 | 2014 | Vert | Microbe | n | Parasite | Endo | horizontal | environmental | resident | 14 | 14 | 11 | 38 | 1.00 | 1.00 | 38 | 1 | 999 | 0.00100 | Parafit |
| Lei_&_Olival_2014 | 2014 | Vert | Microbe | n | Parasite | Endo | horizontal | environmental | resident | 4 | 4 | 4 | 20 | 1.00 | 1.00 | 20 | 1 | 999 | 0.00100 | Parafit |
| Lei_&_Olival_2014 | 2014 | Vert | Microbe | n | Parasite | Endo | horizontal | environmental | resident | 14 | 14 | 12 | 19 | 1.00 | 1.00 | 19 | 1 | 999 | 0.85800 | Parafit |
| Lei_&_Olival_2014 | 2014 | Vert | Microbe | n | Parasite | Endo | horizontal | environmental | resident | 9 | 9 | 8 | 13 | 1.00 | 1.00 | 13 | 1 | 999 | 0.02900 | Parafit |
| Lei_&_Olival_2014 | 2014 | Vert | Microbe | n | Parasite | Endo | horizontal | environmental | resident | 35 | 35 | 22 | 119 | 1.09 | 1.17 | 109 | 1 | 999 | 0.00010 | Parafit |
| Lei_&_Olival_2014 | 2014 | Vert | Microbe | n | Parasite | Endo | horizontal | environmental | resident | 6 | 6 | 5 | 7 | 1.00 | 1.00 | 7 | 1 | 999 | 0.75870 | Parafit |
| LewisRogers_&_Crandall_2009 | 2009 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 6 | 6 | NA | 27 | 1.00 | 1.00 | 27 | 11 | 999 | 0.47000 | Parafit |
| Li_et_al_2017 | 2017 | Plant | Microbe | y | Parasite | Endo | horizontal | vector | Resident | 12 | 12 | 9 | 12 | 1.71 | 2.00 | 7 | 1 | 999 | 0.50505 | Parafit |
| Li_et_al_2017 | 2017 | Plant | Microbe | y | Parasite | Endo | horizontal | autonomous | resident | 12 | 12 | 9 | 12 | 1.71 | 1.57 | 7 | 1 | 999 | 0.50505 | Parafit |
| Li_et_al_2018 | 2018 | Vert | Invert | y | Parasite | Endo | horizontal | trophic | resident | 68 | 68 | 34 | 129 | 1.45 | NA | 89 | 89 | 100000 | 0.00100 | Parafit |
| Light_&_Hafner_2008 | 2008 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 44 | 21 | 4 | 21 | 1.00 | 1.00 | 21 | 1 | 999 | 0.00100 | Parafit |
| Liu_et_al_2016 | 2016 | Plant | Microbe | y | Parasite | Endo | NA | environmental | resident | 13 | 13 | 10 | 44 | 1.57 | 2.71 | 28 | 28 | 9999 | 0.02510 | Parafit |
| Liu_et_al_2016 | 2016 | Plant | Microbe | y | Parasite | Endo | NA | environmental | resident | 19 | 19 | 16 | 76 | 3.30 | 3.22 | 23 | 23 | 9999 | 0.02030 | Parafit |
| Maneesakorn_et_al_2011 | 2011 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 12 | 12 | 1 | 12 | 1.00 | 1.00 | 12 | 1 | 999 | 0.00100 | Parafit |
| Martinez_et_al_2011 | 2011 | Invert | Microbe | y | Parasite | Endo | horizontal | vector | resident | 16 | 16 | 1 | 22 | 2.44 | 1.56 | 9 | 1 | 999 | 0.07500 | Parafit |
| Martinez_San_udo_&_Girolami_2009 | 2009 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 19 | 19 | 10 | 19 | 1.12 | 1.18 | 17 | 5 | 999 | 0.00300 | Parafit |
| Matthews_et_al_2018 | 2018 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 12 | 12 | 6 | 33 | 1.00 | 0.58 | 33 | 1 | 99900 | 0.00400 | Parafit |
| Mattiucci_&_nascetti_2008 | 2008 | Vert | Invert | y | Parasite | Endo | horizontal | trophic | resident | 7 | 7 | 15 | 12 | 1.33 | 1.78 | 9 | 1 | 100 | 0.05000 | Parafit |
| Mazzon_et_al_2010 | 2010 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 17 | 17 | 10 | 17 | 1.06 | 1.13 | 16 | 3 | 999 | 0.00700 | Parafit |
| McFrederick_&_Taylor_2013 | 2013 | Invert | Invert | y | NA | Ecto | vertical | vertical | resident | 7 | 7 | 3 | 7 | 1.00 | 1.00 | 7 | 1 | 10000 | 0.00600 | Parafit |
| McKee_et_al_2016 | 2016 | Vert | Microbe | n | Parasite | Endo | horizontal | vector | resident | 66 | 66 | 42 | 184 | 1.06 | 1.09 | 173 | 1 | 10000 | 0.00010 | Parafit |
| McLeish_&_Van_Noort_2012 | 2012 | Plant | Invert | y | Mutualist | Endo | horizontal | autonomous | resident | 26 | 26 | 1 | 65 | 1.00 | 1.00 | 65 | 6 | 9999 | 0.18000 | Parafit |
| Megia-palma_et_al_2018 | 2018 | Vert | Microbe | y | Parasite | Endo | horizontal | vector | resident | 16 | 16 | 8 | 23 | 1.53 | 1.27 | 15 | 8 | 999 | 0.00400 | Parafit |
| Mehdiabadi_et_al_2012 | 2012 | Invert | Microbe | y | Mutualist | Ecto | both | environmental | resident | 99 | 11 | 1 | 11 | 1.00 | 1.00 | 11 | 1 | 999 | 0.00100 | Parafit |
| Mendlova_et_al_2012 | 2012 | Vert | Invert | y | Parasite | Ecto | horizontal | autonomous | resident | 6 | 6 | 5 | 34 | 1.17 | 1.21 | 29 | 2 | 999 | 0.02400 | Parafit |
| Merville_et_al._2013 | 2013 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 9 | 9 | 1 | 9 | 1.00 | 1.00 | 9 | 1 | 999 | 0.00350 | Parafit |
| Millanes_et_al_2014 | 2014 | Microbe | Microbe | y | Parasite | Endo | both | environmental | resident | 16 | 16 | 2 | 16 | 1.00 | 1.00 | 16 | 1 | 999 | 0.23600 | Parafit |
| Miyaki_et_al_2016 | 2016 | Vert | Microbe | n | Mutualist | Endo | horizontal | environmental | resident | 8 | 8 | 4 | 54 | 3.18 | 2.24 | 17 | 1 | 999 | 0.00100 | Parafit |
| Mondo_et_al_2012 | 2012 | Microbe | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 55 | 5 | 2 | 5 | 1.00 | 1.00 | 5 | 1 | 10000 | 0.00010 | Parafit |
| Nouioui_et_al_2014 | 2014 | Plant | Microbe | n | Mutualist | Endo | horizontal | vector | resident | 9 | 9 | 1 | 20 | 1.00 | 1.00 | 20 | 1 | 9999 | 0.33000 | Parafit |
| Patra_et_al_2018 | 2018 | Vert | Invert | y | Parasite | Endo | horizontal | environmental | resident | 24 | 24 | 24 | 31 | 1.00 | 1.00 | 31 | 1 | 999 | 0.00100 | Parafit |
| Pellissier_et_al_2013 | 2013 | Plant | Invert | y | Parasite | Ecto | horizontal | autonomous | visitor | 104 | 104 | NA | NA | NA | NA | 97 | NA | 10000 | 0.00010 | Parafit |
| Perkins_2010 | 2010 | Vert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 61 | 61 | 50 | 75 | 1.00 | 1.00 | 75 | NA | 9999 | 0.92600 | Parafit |
| Peterson_et_al_2010 | 2010 | Plant | Microbe | y | Parasite | Endo | horizontal | environmental | resident | 11 | 11 | 1 | 25 | 2.08 | 1.58 | 12 | 1 | 9999 | 0.00010 | Parafit |
| Polme_et_al_2014.pdf | 2014 | Plant | Microbe | n | Mutualist | Endo | horizontal | vector | resident | 22 | 22 | 1 | NA | NA | NA | 43 | 1 | 999 | 0.00100 | Parafit |
| Quek_et_al_2004 | 2004 | Plant | Invert | y | Mutualist | Ecto | horizontal | autonomous | resident | 9 | 9 | 1 | 19 | 1.90 | 1.50 | 10 | 1 | 999 | 0.77900 | Parafit |
| Ramasindrazana_et_al_2017 | 2017 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 15 | 15 | 6 | 26 | 2.89 | 1.78 | 9 | 5 | 999 | 0.00100 | Parafit |
| Ricklefs_et_al_2004 | 2004 | Vert | Microbe | y | Parasite | Endo | horizontal | vector | resident | 44 | 44 | 1 | 121 | 1.86 | NA | 65 | 2 | 100 | 0.63000 | Parafit |
| Riess_et_al_2018 | 2018 | Plant | Microbe | Y | Parasite | Endo | horizontal | autonomous | resident | 11 | 11 | 5 | 11 | 1.00 | 1.00 | 11 | 2 | 9999 | 0.26400 | Parafit |
| Savio_2011 | 2011 | Invert | Microbe | n | Mutualist | Endo | horizontal | NA | resident | 17 | 17 | 10 | 17 | 1.00 | 1.00 | 17 | 3 | 999 | 0.00700 | Parafit |
| Schardl_et_al_2008 | 2008 | Plant | Microbe | y | Mutualist | Endo | both | NA | resident | 25 | 25 | 16 | 25 | 0.96 | 1.00 | 26 | 2 | 1000 | 0.00100 | Parafit |
| Sibbald_et_al_2017 | 2017 | Microbe | Microbe | y | Mutualist | Endo | horizontal | environmental | resident | 7 | 7 | 1 | 7 | 1.00 | 1.00 | 7 | 1 | 9999 | 0.00010 | Parafit |
| Simkova_et_al_2013 | 2013 | Vert | Invert | y | Parasite | Ecto | horizontal | autonomous | resident | 5 | 5 | 1 | 21 | 1.00 | 1.00 | 21 | 1 | 999 | 0.01300 | Parafit |
| Singh_et_al_2016 | 2016 | Microbe | Plant | y | Mutualist | Endo | horizontal | environmental | resident | 23 | 23 | 1 | 28 | 1.40 | 1.20 | 20 | 1 | 9999 | 0.00020 | Parafit |
| Singh_et_al_2017 | 2017 | Microbe | Plant | y | Mutualist | Endo | vertical | vertical | resident | 23 | 17 | 1 | 25 | 1.25 | 1.15 | 20 | 1 | 9999 | 0.00020 | Parafit |
| Sontowski_et_al_2015 | 2015 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 14 | 14 | 4 | 14 | 1.00 | 1.00 | 14 | 1 | 100000 | 0.00300 | Parafit |
| Sorenson_et_al_2004 | 2004 | Vert | Vert | y | Parasite | Ecto | horizontal | autonomous | resident | 33 | 33 | 10 | 34 | 1.62 | 1.43 | 21 | 1 | 1000 | 0.00100 | Parafit |
| Souza_et_al_2018 | 2018 | Vert | Microbe | n | Parasite | Endo | both | body fluid | resident | 8 | 8 | 8 | 19 | 1.06 | 1.00 | 18 | 1 | 999 | 0.00500 | Parafit |
| Stireman_et_al_2010 | 2010 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 15 | 15 | 10 | 39 | 1.39 | 1.25 | 28 | 1 | 9999 | 0.00080 | Parafit |
| Sudakaran_et_al_2015 | 2015 | Invert | Microbe | n | Mutualist | Endo | NA | NA | resident | 16 | 16 | 7 | 22 | 1.00 | 1.00 | 22 | 1 | 1000 | 0.97800 | Parafit |
| Sudakaran_et_al_2015 | 2015 | Invert | Microbe | n | Mutualist | Endo | NA | NA | resident | 14 | 14 | 6 | 26 | 1.00 | 1.00 | 26 | 1 | 1000 | 0.97400 | Parafit |
| Summers_&_Rouse_2014 | 2014 | Invert | Invert | y | Parasite | Endo/Ecto | horizontal | autonomous | resident | 53 | 53 | 36 | 78 | 1.13 | 1.23 | 69 | 10 | 9999 | 0.00050 | Parafit |
| Susoy_&_Herrmann_2014 | 2014 | Invert | Invert | y | Mutualist | Ecto | vertical | vertical | resident | 35 | 35 | 7 | 37 | 1.42 | 1.31 | 26 | 1 | 9999 | 0.00010 | Parafit |
| Swafford_&_Bond_2010 | 2010 | Invert | Invert | y | Parasite | Ecto | horizontal | contact | resident | 7 | 7 | 1 | 7 | 1.00 | 1.00 | 7 | 1 | 9999 | 0.31900 | Parafit |
| Sweet_&_Johnson_2018 | 2018 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 13 | 13 | 4 | 14 | 2.80 | 1.80 | 5 | 1 | 100000 | 0.00500 | Parafit |
| Sweet_et_al_2016 | 2016 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 52 | 52 | 25 | 57 | 1.33 | 1.21 | 43 | 1 | 100000 | 0.00001 | Parafit |
| Sweet_et_al_2016 | 2016 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 52 | 52 | 25 | 58 | 1.18 | 1.14 | 49 | 4 | 100000 | 0.00001 | Parafit |
| Sweet_et_al_2017 | 2017 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 12 | 12 | 5 | 15 | 1.15 | 1.23 | 13 | 3 | 999 | 0.06900 | Parafit |
| Sweet_et_al_2018a | 2018 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 259 | 259 | 163 | 283 | 1.63 | 1.61 | 174 | 11 | 9999 | 0.00010 | Parafit |
| Sweet_et_al_2018b | 2018 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 11 | 11 | 4 | 13 | 1.86 | 1.43 | 7 | 1 | 100000 | 0.00500 | Parafit |
| Tao_et_al_2013 | 2013 | Vert | Microbe | n | Parasite | Endo | horizontal | NA | resident | 7 | 7 | 7 | 10 | 1.00 | 1.00 | 10 | 1 | 99999 | 0.17592 | Parafit |
| Toju_et_al_2013 | 2013 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 27 | 27 | 4 | 27 | 1.00 | 1.00 | 27 | 1 | 99999 | 0.00001 | Parafit |
| Urban_&_Cryan_2012 | 2012 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 30 | 30 | 29 | 30 | 1.00 | 1.00 | 30 | 1 | 1000 | 0.00100 | Parafit |
| Urban_&_Cryan_2012 | 2012 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 40 | 40 | 38 | 40 | 1.00 | 1.00 | 40 | 1 | 1000 | 0.00300 | Parafit |
| Viale_et_al_2015 | 2015 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 23 | 23 | 1 | 23 | 1.00 | 1.00 | 23 | 1 | 9999 | 0.00100 | Parafit |
| Won_et_al_2008 | 2008 | Invert | Microbe | n | Mutualist | Endo | horizontal | contact | resident | 15 | 15 | 5 | 15 | 1.00 | 1.00 | 15 | 1 | 100 | 0.37000 | Parafit |
| Xu_et_el_2017 | 2017 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 20 | 20 | 8 | 20 | 1.00 | 1.00 | 20 | 1 | 999 | 0.00100 | Parafit |
| Zhang_et_al_2017 | 2017 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 44 | 44 | 24 | 44 | 1.00 | 1.00 | 44 | 1 | 100000 | 0.00001 | Parafit |
| Badets_et_al_2011 | 2011 | Vert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 17 | 17 | 13 | 17 | 1.00 | 1.00 | 17 | 4 | 1000 | 0.01000 | TreeMap |
| Banks_et_al_2006 | 2006 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 18 | 18 | 6 | 30 | 2.00 | 1.60 | 15 | 2 | 1000 | 0.01000 | TreeMap |
| Bochkov_et_al_2011 | 2011 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 6 | 6 | NA | 9 | 1.00 | 1.00 | 9 | 6 | 100 | 0.01000 | TreeMap |
| Charleston_&_Perkins_2003 | 2003 | Vert | Microbe | y | Parasite | Endo | horizontal | vector | resident | 9 | 9 | 1 | 9 | 1.00 | 1.00 | 9 | 1 | 100 | 0.03000 | TreeMap |
| Charleston_&_Robertson_2002 | 2002 | Vert | Microbe | n | Parasite | Endo | horizontal | body fluid | resident | 12 | 12 | 12 | 12 | 1.00 | 1.00 | 12 | 1 | 1000 | 0.01500 | TreeMap |
| Chauvatcharin_et_al_2006 | 2006 | Microbe | Microbe | n | Parasite | Endo | both | NA | resident | 14 | 14 | 7 | 15 | 1.07 | 1.00 | 14 | 1 | 1000 | 0.33100 | TreeMap |
| Clark_et_al_2000 | 2000 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 17 | 17 | 4 | 17 | 1.00 | 1.00 | 17 | 1 | 1000 | 0.00100 | TreeMap |
| Clayton_&_Johnson_2003 | 2003 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 13 | 13 | 6 | 14 | 1.08 | 1.15 | 13 | 1 | 10000 | 0.00060 | TreeMap |
| Clayton_&_Johnson_2003 | 2003 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 13 | 13 | 6 | 16 | 1.60 | 1.50 | 10 | 1 | 10000 | 0.15300 | TreeMap |
| Cui_et_al_2012 | 2012 | Vert | Microbe | n | Parasite | Endo | horizontal | body fluid | resident | 7 | 7 | 3 | 7 | 1.00 | 1.00 | 7 | 7 | 10000 | 0.36600 | TreeMap |
| Dabert_2001__ | 2001 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 21 | 21 | 12 | 22 | 1.00 | 1.00 | 22 | 9 | 10000 | 0.00100 | TreeMap |
| Deng_et_al_2013 | 2013 | Invert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 7 | 7 | 4 | 10 | 1.00 | 1.00 | 10 | 2 | 999 | 0.53000 | TreeMap |
| Desai_et_al_2010 | 2010 | Microbe | Microbe | n | Mutualist | Ecto | both | contact | resident | 9 | 9 | 2 | 8 | 0.89 | 1.00 | 9 | 1 | 1000 | 0.00100 | TreeMap |
| Desdevises_et_al_2002 | 2002 | Vert | Invert | y | Parasite | Ecto | horizontal | autonomous | resident | 14 | 14 | 11 | 39 | 1.95 | 1.65 | 20 | 2 | 999 | 0.31700 | TreeMap |
| Downie_&_Gullan_2005 | 2005 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 21 | 21 | 16 | 21 | 1.00 | 1.00 | 21 | 1 | 1000 | 0.00100 | TreeMap |
| Erpenbeck_et_al_2002 | 2002 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 6 | 6 | 5 | 6 | 1.00 | 1.00 | 6 | 1 | 10000 | 0.02000 | TreeMap |
| Etherington_et_al_2006 | 2006 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 7 | 7 | 7 | 8 | 1.00 | 1.00 | 8 | 1 | 1000 | 0.01000 | TreeMap |
| Farrell_1998 | 1998 | Plant | Invert | y | Parasite | Ecto | horizontal | autonomous | resident | 21 | 21 | 2 | 21 | 1.00 | 1.00 | 21 | 6 | 1000 | 0.07000 | TreeMap |
| Gottschling_et_al_2011 | 2011 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 43 | 43 | 40 | 78 | 1.00 | 1.00 | 78 | 30 | 1000 | 0.00100 | TreeMap |
| Hendricks_et_al_2013 | 2013 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 19 | 19 | 17 | 20 | 1.25 | 1.25 | 16 | 1 | 1000 | 0.02100 | TreeMap |
| Hosokawa_et_al_2006 | 2006 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 7 | 7 | 3 | 7 | 1.00 | 1.00 | 7 | 1 | 1000 | 0.00100 | TreeMap |
| Hugot_1999 | 1999 | Vert | Invert | y | Parasite | Endo | horizontal | trophic | resident | 10 | 10 | 9 | 11 | 1.00 | 1.00 | 11 | 6 | 1000 | 0.00100 | TreeMap |
| Hugot_et_al_2003 | 2003 | Vert | Microbe | y | Parasite | Endo | horizontal | contact | resident | 19 | 19 | 12 | 19 | 1.00 | 1.00 | 19 | 1 | 1000 | 0.00100 | TreeMap |
| Huyse_&_Volckaert_2005 | 2005 | Vert | Invert | y | Parasite | Ecto | horizontal | autonomous | resident | 8 | 8 | 3 | 22 | 1.29 | 1.29 | 17 | 1 | 100 | 0.02000 | TreeMap |
| IkedaOhtsubo_&_Brune_2009 | 2009 | Microbe | Microbe | n | Mutualist | Ecto | vertical | vertical | resident | 11 | 11 | 1 | 11 | 1.00 | 1.00 | 11 | 1 | 1000 | 0.00100 | TreeMap |
| Jackson_&_Charleston_1994 | 1994 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 10 | 10 | 8 | 10 | 1.00 | 1.00 | 10 | 1 | 100 | 0.19000 | TreeMap |
| Jackson_&_Charleston_1994 | 1994 | Vert | Microbe | n | Parasite | Endo | horizontal | bodily fluid | resident | 10 | 10 | 7 | 10 | 1.00 | 1.00 | 10 | 1 | 100 | 0.01000 | TreeMap |
| Jackson_&_Charleston_1994 | 1994 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 11 | 11 | 7 | 12 | 1.00 | 1.00 | 12 | 1 | 100 | 0.05000 | TreeMap |
| Jackson_&_Charleston_1994 | 1994 | Vert | Microbe | n | Parasite | Endo | horizontal | bodily fluid | resident | 13 | 13 | 9 | 15 | 1.00 | 1.00 | 15 | 1 | 100 | 0.18000 | TreeMap |
| Jackson_&_Charleston_1994 | 1994 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 14 | 14 | 7 | 17 | 1.00 | 1.00 | 17 | 1 | 100 | 0.01000 | TreeMap |
| Jeong_et_al_1999 | 1999 | Plant | Microbe | n | Mutualist | Endo | horizontal | vector | resident | 12 | 12 | 12 | 19 | 1.06 | 1.00 | 18 | 1 | 1000 | 0.23000 | TreeMap |
| Johnson_et_al_2002 | 2002 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 25 | 25 | 23 | 25 | 1.32 | 1.42 | 19 | 19 | 1000 | 0.23000 | TreeMap |
| Johnson_et_al_2003 | 2003 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 28 | 28 | 15 | 31 | 1.48 | 1.33 | 21 | 1 | 100 | 0.03000 | TreeMap |
| Johnson_et_al_2006 | 2006 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 10 | 10 | NA | 11 | 1.00 | 1.00 | 11 | 11 | 1000 | 0.40200 | TreeMap |
| Jousselin_et_al_2008 | 2009 | Plant | Invert | y | Mutualist | Endo | horizontal | autonomous | resident | 15 | 15 | 1 | 15 | 1.07 | 1.07 | 14 | 6 | 10000 | 0.01000 | TreeMap |
| Jousselin_et_al_2008 | 2009 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 13 | 13 | 1 | 14 | 1.00 | 1.00 | 14 | 1 | 10000 | 0.01000 | TreeMap |
| Jousselin_et_al_2008 | 2009 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 13 | 13 | 1 | 13 | 1.00 | 1.00 | 13 | 1 | 10000 | 0.01000 | TreeMap |
| Jousselin_et_al_2008 | 2009 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 16 | 16 | 1 | 18 | 1.06 | 1.06 | 17 | 2 | 10000 | 0.01000 | TreeMap |
| Jousselin_et_al_2009 | 2009 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 55 | 22 | 1 | 22 | 1.00 | 1.00 | 22 | 1 | 10000 | 0.00100 | TreeMap |
| Kawaida_et_al_2013 | 2013 | Invert | Plant | y | Mutualist | Endo | vertical | vertical | resident | 6 | 6 | 1 | 6 | 1.00 | 1.00 | 6 | 1 | 10000 | 0.00350 | TreeMap |
| Kawakita_et_al_2004 | 2004 | Plant | Invert | y | Mutualist | Ecto | horizontal | autonomous | visitor | 18 | 18 | 1 | 18 | 1.00 | 1.00 | 18 | 1 | 999 | 0.01900 | TreeMap |
| Kelley_&_Farrell_1998 | 1998 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 41 | 41 | 1 | 89 | 6.85 | 1.92 | 13 | 1 | 100 | 0.28000 | TreeMap |
| Kikuchi_et_al_2009 | 2009 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 14 | 14 | 5 | 14 | 1.00 | 1.00 | 14 | 1 | 1000 | 0.00100 | TreeMap |
| Lanterbecq_et_al_2010 | 2010 | Invert | Invert | y | Parasite | Endo/Ecto | horizontal | contact | resident | 16 | 16 | 12 | 16 | 1.00 | 1.00 | 16 | 5 | 5000 | 0.04000 | TreeMap |
| Light_&_Hafner_2008 | 2008 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 44 | 21 | 4 | 21 | 1.00 | 1.00 | 21 | 1 | 1000 | 0.00100 | TreeMap |
| LimFong_et_al_2008 | 2008 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 5 | 5 | 1 | 5 | 1.00 | 1.00 | 5 | 2 | 1000 | 0.11000 | TreeMap |
| Liu_et_al_2013 | 2013 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 37 | 37 | 1 | 37 | 1.00 | 1.00 | 37 | 1 | 1000 | 0.00100 | TreeMap |
| Liu_et_al_2014 | 2014 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 27 | 27 | 3 | 29 | 1.00 | 1.00 | 29 | 1 | 1000 | 0.01000 | TreeMap |
| LopezVaamonde_et_al_2001 | 2001 | Invert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 15 | 15 | 1 | 15 | 1.00 | 1.00 | 15 | 1 | 1000 | 0.00100 | TreeMap |
| LopezVaamonde_et_al_2003 | 2003 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 33 | 33 | 33 | 77 | 1.00 | 1.00 | 77 | 1 | 1000 | 0.21300 | TreeMap |
| LopezVaamonde_et_al_2005 | 2005 | Invert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 28 | 28 | 1 | 35 | 2.33 | 1.33 | 15 | 1 | 1000 | 0.24800 | TreeMap |
| Martin_et_al_1999 | 1999 | Vert | Microbe | n | Parasite | Endo | horizontal | body fluid | resident | 38 | 38 | 38 | 48 | 1.00 | 1.00 | 48 | 1 | 1000 | 0.00100 | TreeMap |
| Martin_et_al_2003 | 2003 | Vert | Microbe | n | Parasite | Endo | horizontal | body fluid | resident | 14 | 14 | 14 | 16 | 1.00 | 1.00 | 16 | 1 | 100 | 0.01000 | TreeMap |
| Martin_et_al_2003 | 2003 | Vert | Microbe | n | Parasite | Endo | horizontal | body fluid | resident | 9 | 9 | 9 | 9 | 1.00 | 1.00 | 9 | 1 | 100 | 0.05000 | TreeMap |
| Martin_et_al_2003 | 2003 | Vert | Microbe | n | Parasite | Endo | horizontal | body fluid | resident | 17 | 17 | 17 | 23 | 1.00 | 1.00 | 23 | 1 | 100 | 0.21000 | TreeMap |
| Mazzon_et_al_2010 | 2010 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 17 | 17 | 10 | 17 | 1.06 | 1.13 | 16 | 3 | 1000 | 0.00100 | TreeMap |
| Morelli_&_Spicer_2007 | 2007 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 6 | 6 | 6 | 6 | 1.00 | 1.00 | 6 | 1 | 10000 | 0.00950 | TreeMap |
| Muniz_et_al_2013 | 2013 | Vert | Microbe | n | Parasite | Endo | horizontal | body fluid | resident | 17 | 17 | 11 | 17 | 1.00 | 1.00 | 17 | 1 | 10000 | 0.00001 | TreeMap |
| Musser_et_al_2010 | 2010 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 6 | 6 | 1 | 6 | 1.00 | 1.00 | 6 | 6 | 100 | 0.30000 | TreeMap |
| Pagan_et_al_2010 | 2010 | Plant | Microbe | n | Parasite | Endo | horizontal | contact | resident | 10 | 10 | 10 | 13 | 1.00 | 1.00 | 13 | 1 | 1000 | 0.24000 | TreeMap |
| Page_1996 | 1996 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 15 | 15 | 6 | 17 | 1.00 | 1.00 | 17 | 2 | 1000 | 0.00100 | TreeMap |
| Page_et_al_1998 | 1998 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 7 | 7 | 1 | 8 | 1.00 | 1.00 | 8 | 1 | 100 | 0.01000 | TreeMap |
| Page_et_al_2004 | 2004 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 12 | 12 | 4 | 12 | 1.00 | 1.00 | 12 | 1 | 1000 | 0.00100 | TreeMap |
| Page_et_al_2004 | 2004 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 11 | 11 | 5 | 13 | 1.00 | 1.00 | 13 | 4 | 100 | 0.25000 | TreeMap |
| Page_et_al_2004 | 2004 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 13 | 13 | 5 | 14 | 1.00 | 1.00 | 14 | 1 | 100 | 0.46000 | TreeMap |
| Page_et_al_2004 | 2004 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 9 | 9 | 3 | 9 | 1.00 | 1.00 | 9 | 1 | 100 | 0.36000 | TreeMap |
| Paterson_&_Poulin_1999 | 1999 | Vert | Invert | y | Parasite | Ecto | horizontal | autonomous | resident | 8 | 8 | 8 | 12 | 1.20 | 2.20 | 10 | 1 | 100 | 0.01000 | TreeMap |
| Paterson_et_al_2000 | 2000 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 11 | 11 | 5 | 14 | 1.00 | 1.00 | 14 | 5 | 100 | 0.01000 | TreeMap |
| Percy_et_al_2004 | 2004 | Plant | Invert | y | Parasite | Ecto | both | autonomous | resident | 35 | 35 | 8 | 56 | 1.22 | 1.09 | 46 | 4 | 1000 | 0.00500 | TreeMap |
| PerezLosada_et_al_2006 | 2006 | Vert | Microbe | n | Parasite | Endo | horizontal | NA | resident | 9 | 9 | 9 | 11 | 1.00 | 1.00 | 11 | 1 | 100 | 0.01000 | TreeMap |
| Perlman_et_al_2003 | 2003 | Invert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 16 | 16 | 1 | 17 | 1.89 | 1.33 | 9 | 1 | 1000 | 0.14000 | TreeMap |
| Quek_et_al_2004 | 2004 | Plant | Invert | y | Mutualist | Ecto | horizontal | autonomous | resident | 9 | 9 | 1 | 9 | 1.00 | 1.00 | 6 | 1 | 1000 | 1.00000 | TreeMap |
| Ramsden_et_al_2008 | 2008 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 33 | 33 | 20 | 38 | 1.00 | 1.00 | 38 | 1 | 1000 | 1.00000 | TreeMap |
| Reed_et_al_2007 | 2007 | Vert | Invert | y | Parasite | Ecto | horizontal | contact | resident | 4 | 4 | 4 | 6 | 1.20 | 1.00 | 5 | 4 | 1000 | 0.05000 | TreeMap |
| Refregier_et_al_2008 | 2008 | Plant | Microbe | y | Parasite | Endo | horizontal | vector | resident | 18 | 18 | 7 | 20 | 1.25 | 1.25 | 16 | 1 | 3000 | 0.50000 | TreeMap |
| Riess_et_al_2018 | 2018 | Plant | Microbe | Y | Parasite | Endo | horizontal | autonomous | resident | 11 | 11 | 5 | 11 | 1.00 | 1.00 | 11 | 2 | 1000 | 0.33500 | TreeMap |
| Shoemaker_et_al_2002 | 2002 | Invert | Microbe | n | Mutualist | Endo | both | NA | resident | 20 | 20 | 7 | 23 | 1.00 | 1.00 | 23 | 1 | 100 | 0.47000 | TreeMap |
| Simkova_et_al_2013 | 2013 | Vert | Invert | y | Parasite | Ecto | horizontal | autonomous | resident | 5 | 5 | 1 | 21 | 1.00 | 1.00 | 21 | 1 | 999 | 0.04500 | TreeMap |
| Six_&_Paine_1999 | 1999 | Microbe | Invert | y | Mutualist | Ecto | horizontal | vector | resident | 6 | 6 | 1 | 6 | 1.00 | 1.00 | 6 | 3 | 1000 | 0.03100 | TreeMap |
| Skerikova_et_al_2001 | 2001 | Vert | Invert | y | Parasite | Endo | horizontal | trophic | resident | 7 | 7 | 7 | 7 | 1.00 | 1.00 | 7 | 1 | 10000 | 0.40000 | TreeMap |
| Smith_et_al_2008b | 2008 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 20 | 20 | 3 | 20 | 1.00 | 1.00 | 20 | 14 | 100 | 0.05000 | TreeMap |
| Sorenson_et_al_2004 | 2004 | Vert | Vert | y | Parasite | Ecto | horizontal | autonomous | resident | 33 | 33 | 10 | 34 | 1.62 | 1.43 | 21 | 1 | 100 | 0.15000 | TreeMap |
| Subbotin_et_al_2004 | 2004 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 16 | 16 | 16 | 21 | 1.00 | 1.00 | 21 | 4 | 1000 | 0.00100 | TreeMap |
| Switzer_et_al_2005 | 2005 | Vert | Microbe | n | Parasite | Endo | horizontal | body fluid | resident | 46 | 46 | 17 | 51 | 1.00 | 1.00 | 51 | 1 | 10000 | 0.00700 | TreeMap |
| Urban_&_Cryan_2012 | 2012 | Invert | Microbe | n | Mutualist | Endo | vertical | autonomous | resident | 40 | 40 | 38 | 40 | 1.00 | 1.00 | 40 | 1 | 1000 | 0.00100 | TreeMap |
| Urban_&_Cryan_2012 | 2012 | Invert | Microbe | n | Mutualist | Endo | vertical | autonomous | resident | 30 | 30 | 29 | 30 | 1.00 | 1.00 | 30 | 1 | 1000 | 0.00100 | TreeMap |
| Vanhove_et_al_2015 | 2015 | Vert | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 19 | 19 | 10 | 28 | 1.00 | 1.00 | 28 | 1 | 10000 | 0.04210 | TreeMap |
| Weckstein_2004 | 2004 | Vert | Invert | y | Parasite | Ecto | both | contact | resident | 11 | 11 | 1 | 11 | 2.20 | 1.80 | 5 | 1 | 10000 | 0.89000 | TreeMap |
| Weiblen_&_Bush_2002 | 2002 | Plant | Invert | y | Mutualist | Endo | horizontal | autonomous | resident | 19 | 19 | 1 | 19 | 1.00 | 1.00 | 19 | 1 | 10000 | 0.01950 | TreeMap |
| Weiblen_&_Bush_2002 | 2002 | Plant | Invert | y | Parasite | Endo | horizontal | autonomous | resident | 12 | 12 | 1 | 18 | 1.00 | 1.00 | 18 | 1 | 10000 | 0.12150 | TreeMap |
| Wu_et_al_2008 | 2008 | Plant | Microbe | n | Parasite | Endo | horizontal | vector | resident | 8 | 8 | 8 | 10 | 1.00 | 1.00 | 10 | 1 | 100 | 0.01000 | TreeMap |
| Xu_et_al_2017 | 2017 | Invert | Microbe | n | Mutualist | Endo | vertical | vertical | resident | 20 | 20 | 8 | 20 | 1.00 | 1.00 | 20 | 1 | 1000 | 0.00100 | TreeMap |
| Yan_et_al_2011 | 2011 | Vert | Microbe | n | Parasite | Endo | horizontal | contact | resident | 8 | 8 | 8 | 15 | 1.00 | 1.00 | 15 | 1 | 999 | 0.98900 | TreeMap |
A. authors: The authors of the study and the date (citation form).
B. year: The year of publication of the study.
C. host_tax_broad: Separation of the host group according to broader taxonomic units (e.g. vertebrate, invertebrate, microbe, plant).
D. symbiont_tax_broad: Separation of the symbiont group according to broader taxonomic units (e.g. vertebrate, invertebrate, microbe, plant).
E. symbiont_euk*: Whether the symbiont is eukaryotic (state =‘yes’), or prokaryotic (state=‘no’).
F. symbiosis: The type of symbiont (e.g. parasite or mutualist). For this we followed the definition used by the authors of the study.
G. endo_or_ecto: Whether the symbiont lives outside the host (i.e. is an ectosymbiont), or inside the host (i.e. is an endosymbiont).
H. mode_of_transmission_broad: Whether the symbiont is transmitted vertically, horizontally, or both. For this, we followed the route of transmission specified by the authors of the study.
I. mode_of_transmission_fine*: A finer-scale description of the mode of transmission of the symbiont (e.g. contact, vector, bodily fluid, vertical, trophic).
J. Visiting_symbiont?* Whether the symbiont is resident on the host (resident), or makes visits to the host or hosts (visitor).
K. host_tips_linked: The number of individual host taxa included in the cophylogenetic analysis.
L. host_tips_linked_corrected The same measure as for column K, ‘host_tips_linked,’ but reduced to only include one member of each host species. This is included because some authors include multiple individuals of the same host species. Without correction, this artificially increases the apparent number of host species included in the study.
M. host_genera: A count of the number of host genera included in the cophylogenetic analysis.
N. total_host_symbiont_links The total number of links between host and symbiont taxa recorded in a study. If all symbionts were strict specialists, this would equal the number of symbionts included in the study. However, because symbionts are often associated with more than one host, this value is often higher than the total number of symbionts included in the study.
O. host_range_link_ratio: An estimation of symbiont host specificity, calculated by dividing the total number of links between hosts and symbionts (i.e. ‘total_host_symbiont_links,’ column N), by the total number of symbionts included in the study (i.e. ‘symbiont_tips_linked,’ column Q).
P. host_range_taxonomic_breadth: An alternative estimation of symbiont host specificity, calculated by first summing the number of host taxonomic ranks linked to each symbiont (i.e. single host species = 1, multiple host species in the same genus = 2, multiple host genera = 3, multiple host families = 4, multiple host orders = 5), and dividing by the total number of symbionts included in the study (i.e. ‘symbiont_tips_linked,’ column Q).
Q. symbiont_tips_linked The number of individual symbiont taxa included in the cophylogenetic analysis.
R. symbiont_genera: A count of the number of symbiont genera included in the cophylogenetic analysis.
S. no_randomizations: The number of phylogenetic randomizations performed during the cophylogenetic analysis.
T. p_value: The p-value reported for the cophylogenetic analysis, representing the likelihood that host and symbiont phylogenies display cospeciation.
U. method: Whether TreeMap or ParaFit was used to obtain the reported p value.
Below we present our sample sizes for the two separate methods: TreeMap (Page 1994) and ParaFit (Legendre et al. 2002) (and combined), in terms of effect sizes, papers, and different levels of categorical variables (factors).
# selecting out variables, which we used for our analysis
dat <- full_data %>% select(-symbiont_euk, -mode_of_transmission_fine, -`Visiting_symbiont?`)
pair <- full_pair %>% select(-symbiont_euk, -mode_of_transmission_fine, -`Visiting_symbiont?`)
# making a table of sample sizes for different variables
dat %>% group_by(method) %>%
summarise(
`Effect sizes (analyses)` = n(),
Studies = n(),
Papers = n_distinct(authors),
`Vertebrate hosts` = sum(host_tax_broad == "Vert", na.rm = T), # na.rm is important when NA exists
`Invertebrate hosts` = sum(host_tax_broad == "Invert", na.rm = T),
`Plant hosts` = sum(host_tax_broad == "Plant", na.rm = T),
`Microbe hosts` = sum(host_tax_broad == "Microbe", na.rm = T),
`Vertebrate symbionts` = sum(symbiont_tax_broad == "Vert", na.rm = T),
`Invertebrate symbionts` = sum(symbiont_tax_broad == "Invert", na.rm = T),
`Plant symbionts` = sum(symbiont_tax_broad == "Plant", na.rm = T),
`Microbe symbionts` = sum(symbiont_tax_broad == "Microbe", na.rm = T),
`Parastic relationships` = sum(symbiosis == "Parasite", na.rm = T),
`Mutualistic relatioships` = sum(symbiosis == "Mutualist", na.rm = T),
`Ecto-symbionts` = sum(endo_or_ecto == "Ecto", na.rm = T),
`Endo-symbionts` = sum(endo_or_ecto == "Endo", na.rm = T),
`Ecto/endo-symbionts` = sum(endo_or_ecto == "Endo/Ecto", na.rm = T),
`Horizontal transmission` = sum(mode_of_transmission_broad == "horizontal", na.rm = T),
`Vertical transmission` = sum(mode_of_transmission_broad == "vertical", na.rm = T),
`Horizontal/vertical-transmission` = sum(mode_of_transmission_broad == "both", na.rm = T)
) -> n_table1
# transposing the table and creating that table and adding a correct number of the papers for `Combined`
n_authors <- n_distinct(dat$authors) # the total number of papers
dat$studies <- paste0(dat$authors, dat$host_tax_fine, dat$symbiont_tax_fine, dat$total_host_symbioint_links) # the number of analyses (studies)
n_studies <- n_distinct(dat$studies)
n_table2 <-t(n_table1[,-1])
colnames(n_table2) <- n_table1$method
n_table2 %>% as_tibble(rownames = "Number") %>%
mutate(Combined = Parafit + TreeMap, Combined = replace(Combined, c(2,3), c(n_studies, n_authors))) %>%
rename("Number of" = "Number", "ParaFit (n)" = "Parafit", "TreeMap (n)" = "TreeMap", "Combined (n)" = "Combined") %>%
kable() %>% kable_styling("striped", position = "left") %>%
scroll_box(width = "100%", height = "250px")
| Number of | ParaFit (n) | TreeMap (n) | Combined (n) |
|---|---|---|---|
| Effect sizes (analyses) | 140 | 93 | 233 |
| Studies | 140 | 93 | 211 |
| Papers | 118 | 78 | 180 |
| Vertebrate hosts | 60 | 51 | 111 |
| Invertebrate hosts | 39 | 20 | 59 |
| Plant hosts | 31 | 18 | 49 |
| Microbe hosts | 10 | 4 | 14 |
| Vertebrate symbionts | 1 | 1 | 2 |
| Invertebrate symbionts | 62 | 49 | 111 |
| Plant symbionts | 3 | 1 | 4 |
| Microbe symbionts | 74 | 42 | 116 |
| Parastic relationships | 91 | 70 | 161 |
| Mutualistic relatioships | 48 | 23 | 71 |
| Ecto-symbionts | 41 | 34 | 75 |
| Endo-symbionts | 97 | 58 | 155 |
| Ecto/endo-symbionts | 2 | 1 | 3 |
| Horizontal transmission | 84 | 53 | 137 |
| Vertical transmission | 28 | 15 | 43 |
| Horizontal/vertical-transmission | 23 | 25 | 48 |
#pander(split.cell = 40, split.table = Inf) # not as nice as kable
Note that the numbers of studies and papers does not add up (TreeMap + ParaFit \(\neq\) Combined), because 22 analyses and 16 papers used both the TreeMap and ParaFit methods (the term “papers” here is our variable authors)
Below, we present the number of instances of missing data (cells) for all variables used in our meta-analysis.
# summaring missingness in our dataset
# funs(sum(is.na(.))) needs to be in funs as is.na has "." = each column
dat %>% summarise_all(~sum(is.na(.))) %>% # map(~sum(is.na(.)) # this is an alterantive way
t() %>% as_tibble(rownames = "Variable") %>%
rename("Number of missing data (n)" = "V1") %>%
#pander(split.cell = 40, split.table = Inf)
kable() %>% kable_styling("striped", position = "left") %>%
scroll_box(width = "60%", height = "250px")
| Variable | Number of missing data (n) |
|---|---|
| authors | 0 |
| year | 0 |
| host_tax_broad | 0 |
| symbiont_tax_broad | 0 |
| symbiosis | 1 |
| endo_or_ecto | 0 |
| mode_of_transmission_broad | 5 |
| host_tips_linked | 0 |
| host_tips_linked_corrected | 0 |
| host_genera | 6 |
| total_host_symbioint_links | 3 |
| host_range_link_ratio | 3 |
| host_range_taxonomic_breadth | 7 |
| symbiont_tips_linked | 0 |
| symbiont_genera | 9 |
| no_randomizations | 0 |
| p_value | 0 |
| method | 0 |
| studies | 0 |
# an alternative method using the mi package
#missing_data_tbl <- missing_data.frame(as.data.frame(data))
#show(missing_data_tbl)
We created our effect size (correlation coefficient r and its Fisher’s z transformation Zr) from p values and associated sample sizes (Rosenthal & Rubin 2003). We used half of the“effective sample size” for each effect size (an indicator of congruence):
\[
N_{effective} = \frac{4n_{h}n_{s}} {n_{h} + n_{s}},
\] where \(n_h\) and \(n_s\) are sample size for hosts and symbionts (in our case: host_tips_linked_corrected and symbiont_tips_linked). When \(n = n_h = n_s\), \(N_{effective}\) simplifies to \(2n\). We used \(N_{effective}/2\), which we called the number of “effective pairs” (\(N_{ep}\)) as our sample size.
Also, we created a column with a unique ID for each observation (i.e. an observation level random effect), termed observation, which is required for the rma.mv function in metafor (Viechtbauer 2010).
dat %<>% # getting sample size & observation level random effect
mutate(., sample_size = 2 * (host_tips_linked_corrected * symbiont_tips_linked)/(host_tips_linked_corrected +
symbiont_tips_linked), observation = factor(1:nrow(.))) #, signficance = if_else(p_value <= 0.05, 'y', 'n') )
# making p = 1 to p = (no_randomization - 1)/no_randomization as p = 1 produces t
# value = Inf
dat$p_value <- ifelse(dat$p_value != 1, dat$p_value, (dat$no_randomizations - 1)/dat$no_randomizations)
# calculating effect size
dat %<>% p_to_Zr(p_value, sample_size)
# getting sample size for spp sum(dat$sample_size[match(unique(dat$studies),
# dat$studies)])
# getting effect sizes for pair data
pair %<>% # getting sample size & observation level random effect
mutate(., sample_size = 2 * (host_tips_linked_corrected * symbiont_tips_linked)/(host_tips_linked_corrected +
symbiont_tips_linked), observation = factor(1:nrow(.)))
# making p = 1 to p = (no_randomization - 1)/no_randomization as p = 1 produces t
# value = Inf
pair$p_value_tree <- ifelse(pair$p_value_tree == 1, 0.999, pair$p_value_tree)
# calculating effect size
pair %<>% p_to_Zr(p_value_para, sample_size) %>% rename(rval_para = rval, Zr_para = Zr,
VZr_para = VZr) %>% p_to_Zr(p_value_tree, sample_size) %>% rename(rval_tree = rval,
Zr_tree = Zr, VZr_tree = VZr)
In 16 studies, authors used both ParaFit and TreeMap on the same cophylogeny dataset, but in some of these studies, they used different numbers of randomization between the two methods. Among these 16, 11 studies used ParaFit and TreeMap with the same or similar number of randomization (e.g., 999 and 1,000). We here correlated effect sizes (Zr) between these two methods. In the following analyses, we assume the p-value based effect size, Zr~equivalent (or just Zr) is equivalent for both methods. In such a case, we would require a tight correlation between the method.
# plotting correlations but not displayed A)
cor_1 <- round(with(pair, cor(Zr_para, Zr_tree)), 3)
cor_t1 <- with(pair, cor.test(Zr_para, Zr_tree))
all_plot <- ggplot(pair, aes(Zr_para, Zr_tree)) + geom_point() + geom_smooth(method = "lm") +
labs(x = "Zr (ParaFit)", y = "Zr (TreeMap)") + annotate("text", x = 1, y = -1,
label = paste("r = ", cor_1))
# B)
cor_2 <- round(with(pair[pair$match == "y", ], cor(Zr_para, Zr_tree)), 3)
cor_t2 <- with(pair[pair$match == "y", ], cor.test(Zr_para, Zr_tree))
part_plot <- ggplot(pair[pair$match == "y", ], aes(Zr_para, Zr_tree)) + geom_point() +
geom_smooth(method = "lm") + labs(x = "Zr (ParaFit)", y = "Zr (TreeMap)") + annotate("text",
x = 1, y = -1, label = paste("r = ", cor_2))
comb_plot <- all_plot + part_plot + plot_annotation(tag_levels = "A", tag_suffix = ")")
# comb_plot
There are strong and significant correlations in effect sizes (Zr) between ParaFit and TreeMap (for all the 16 studies, r = 0.737, t = 4.075, df = 14, p = 0.0011, and for the 11 perfectly matching studies, r = 0.781, t = 3.757, df = 9, p = 0.0045). Therefore, we have made an assumption that effect sizes based on p values from these two methods can be equated in our analyses below.
First, we checked what random effects should be put into the main model. To do this we fitted two random effects, authors (i.e. study IDs) and observation; the former term was added to account for non-independence of effect sizes originating from the same papers (i.e., authors).
# 2 random effects & model AIC note that probably only base stuff works outside
# of main chunck so need to create AIC here
ma_test1 <- rma.mv(yi = Zr, V = VZr, random = list(~1 | authors, ~1 | observation),
data = dat)
aic1 <- AIC(ma_test1)
# 1 random effect & model AIC
ma_test2 <- rma.mv(yi = Zr, V = VZr, random = ~1 | authors, data = dat)
aic2 <- AIC(ma_test2)
The model (ma_test1), which included both random factors, had a larger AIC value (279.99) than the model with only one random effect (277.99) . This is because observation hardly accounted for any variance (< 0.0001) compared to authors (0.0811). Therefore, we only included authors as our random factor in subsequent analyses.
We ran intercept models (meta-analyses) with 3 different datasets (ParaFit, TreeMap and both combined; see the explanation of method above). Also, we note that we used adjustments for test statistics and confidence intervals (test = "t"), which is similar to (but not the same as) those proposed by Knapp and Hartung (Knapp & Hartung 2003); probably this approach is more conservative.
# think about making this into a tibble meta-analysis with Parafit
ma_parafit <- rma.mv(yi = Zr, V = VZr, random = ~1 | authors, test = "t", subset = which(method ==
"TreeMap"), data = dat)
# meta-analysis with TreeMap
ma_treemap <- rma.mv(yi = Zr, V = VZr, random = ~1 | authors, test = "t", subset = which(method ==
"Parafit"), data = dat)
# meta-analysis with all the data combined
ma_all <- rma.mv(yi = Zr, V = VZr, test = "t", random = ~1 | authors, data = dat)
Overall effects (meta-analytic means) and 95% confidence intervals (CIs) both in Zr and r, and variance components (V) and heterogeneity, I2 (I2) (Higgins et al. 2003) from the metafor model using the 3 datasets (ParaFit, TreeMap and both combined, or All). Note that in these models, I2[total] = I2[authors] (see (Nakagawa & Santos 2012; Senior et al. 2016)), as we only have one random factor.
# getting I2 for the models could use map()
i2_treemap <- I2(ma_treemap)
i2_parafit <- I2(ma_parafit)
i2_all <- I2(ma_all)
# creating a table
tibble(Dataset = c("Parafit (Zr)", "TreeMap (Zr)", "All (Zr)", "Parafit (r)", "TreeMap (r)",
"All (r)"), `Overall mean` = c(ma_parafit$b, ma_treemap$b, ma_all$b, tanh(ma_parafit$b),
tanh(ma_treemap$b), tanh(ma_all$b)), `Lower CI [0.025]` = c(ma_parafit$ci.lb,
ma_treemap$ci.lb, ma_all$ci.lb, tanh(ma_parafit$ci.lb), tanh(ma_treemap$ci.lb),
tanh(ma_all$ci.lb)), `Upper C [0.975]` = c(ma_parafit$ci.ub, ma_treemap$ci.ub,
ma_all$ci.ub, tanh(ma_parafit$ci.ub), tanh(ma_treemap$ci.ub), tanh(ma_all$ci.ub)),
`V[authors]` = c(ma_parafit$sigma2, ma_treemap$sigma2, ma_all$sigma2, rep(NA,
3)), `I2[total]` = c(i2_parafit[1], i2_treemap[1], i2_all[1], rep(NA, 3))) %>%
kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Dataset | Overall mean | Lower CI [0.025] | Upper C [0.975] | V[authors] | I2[total] |
|---|---|---|---|---|---|
| Parafit (Zr) | 0.518 | 0.425 | 0.610 | 0.084 | 0.526 |
| TreeMap (Zr) | 0.544 | 0.477 | 0.610 | 0.072 | 0.627 |
| All (Zr) | 0.536 | 0.479 | 0.593 | 0.081 | 0.612 |
| Parafit (r) | 0.476 | 0.401 | 0.544 | NA | NA |
| TreeMap (r) | 0.496 | 0.444 | 0.544 | NA | NA |
| All (r) | 0.490 | 0.446 | 0.532 | NA | NA |
These models all gave consistent results including heterogeneity. Given these results, we proceeded with only analyzing the whole dataset (All) from here on.
# https://stackoverflow.com/questions/41919023/ggplot-adding-image-on-top-right-in-two-plots-with-different-scales
# how to add png files to the figure (above) reading image
image_mutualism <- readPNG(here("images/mutualism_transparentbg.png"))
image_parasitism <- readPNG(here("images/parasitism_transparentbg.png"))
# creating a table of results
pred_ma <- get_pred(ma_all)
effect_ma <- get_est(ma_all) %>% left_join(pred_ma)
# creating An orchard plot
fig_ma <- ggplot(data = effect_ma, aes(x = tanh(estimate), y = "Overall mean")) +
scale_x_continuous(limits = c(-1, 1), breaks = seq(-1, 1, by = 0.2)) + geom_quasirandom(data = dat,
aes(x = tanh(Zr), y = "Overall mean", size = (1/VZr) + 3), groupOnX = FALSE,
alpha = 0.2) + # precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F,
size = 0.5, alpha = 0.6) + # CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F,
size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) + # creating dots and different size (bee-swarm and bubbles)
geom_point(size = 3, shape = 21, fill = "black") + annotate("text", x = 0.93, y = 1.15,
label = paste("italic(k)==", length(dat$Zr)), parse = TRUE, hjust = "left", size = 3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n),
" (# of species pair)"))) + theme_bw() + theme(legend.position = c(0, 1),
legend.justification = c(0, 1)) + theme(legend.direction = "horizontal") + # theme(legend.background = element_rect(fill = 'white', colour = 'black')) +
theme(legend.background = element_blank()) + theme(axis.text.y = element_text(size = 10,
colour = "black", hjust = 0.5, angle = 90)) + annotation_custom(rasterGrob(image_mutualism),
xmin = -1.1, xmax = -0.9, ymin = 0.6, ymax = 1.2) + annotation_custom(rasterGrob(image_parasitism),
xmin = -0.9, xmax = -0.7, ymin = 0.6, ymax = 1.2)
# ggsave(plot = fig_ma, filename = 'fig_2a.pdf', height = 2, width = 8) ggploty 0
# does not work (Error in unique.default(x) : unimplemented type 'expression' in
# 'HashTableSetup')
fig_ma
# for Fig 3
a <- ggplot(data = effect_ma, aes(x = tanh(estimate), y = "Overall mean")) + scale_x_continuous(limits = c(-1,
1), breaks = seq(-1, 1, by = 0.2)) + geom_quasirandom(data = dat, aes(x = tanh(Zr),
y = "Overall mean", size = (1/VZr) + 3), groupOnX = FALSE, alpha = 0.2) + # precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F,
size = 0.5, alpha = 0.6) + # CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F,
size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) + # creating dots and different size (bee-swarm and bubbles)
geom_point(size = 3, shape = 21, fill = "black") + annotate("text", x = 0.93, y = 1.15,
label = paste("italic(k)==", length(dat$Zr)), parse = TRUE, hjust = "left", size = 3.5) +
labs(x = "", y = "", size = expression(paste(italic(n), " (# of species pairs)")),
tag = "a") + theme_bw() + theme(legend.position = c(0, 1), legend.justification = c(0,
1)) + theme(legend.direction = "horizontal") + # theme(legend.background = element_rect(fill = 'white', colour = 'black')) +
theme(legend.background = element_blank()) + theme(axis.text.y = element_text(size = 10,
colour = "black", hjust = 0.5, angle = 90)) + annotation_custom(rasterGrob(image_mutualism),
xmin = -1.1, xmax = -0.9, ymin = 0.6, ymax = 1.2) + annotation_custom(rasterGrob(image_parasitism),
xmin = -0.9, xmax = -0.7, ymin = 0.6, ymax = 1.2)
Figure 3a: An orchard plot showing the meta-analytic mean (mean effect size) with its 95% confidence interval (thick line) and 95% prediction interval (thin line), with observed effect sizes based on various sample sizes.
We ran a univariate meta-regression model for each of the following moderators: 1) symbiosis, 2) host_tax_broad, 3) symbiont_tax_broad, 4) host_range_link_ratio, 5) host_range_taxonomic_breadth, 6) mode_of_transmission_broad, and 7) endo_or_ecto. The results from these models are presented in the main text.
In addition to these, we ran three more univariate models: 1) host_tax_symbiosis (equivalent to the interaction term between symbiosis and host_tax_symbiosis; symbiosis*host_tax_symbiosis), 2) symbiont_tax_symbiosis (symbiosis*symbiont_tax_broad), 3) host_symbiont_tax (host_tax_symbiosis*symbiont_tax_broad) and 4) symbiosis_transmission (symbiosis*mode_of_transmission_broad). These moderators are created below:
dat %<>%
# host_tax_broad*symbiosis (host_tax_symbiosis)
mutate(host_tax_symbiosis = str_c(host_tax_broad, symbiosis),
host_tax_symbiosis = ifelse(host_tax_symbiosis == "InvertNA", NA, host_tax_symbiosis),
host_tax_symbiosis = factor(host_tax_symbiosis),
# symbiont_tax_broad*symbiosis (symbiont_tax_symbiosis)
symbiont_tax_symbiosis = factor(str_c(symbiont_tax_broad, symbiosis)),
# host_tax_broad*symbiont_tax_broad (host_symbiont_tax)
host_symbiont_tax = factor(str_c(host_tax_broad, symbiont_tax_broad)),
# symbiosis*mode_of_transmission_broad (symbiosis_transmission)
symbiosis_transmission = factor(str_c(symbiosis, mode_of_transmission_broad)),
# whether p values were the smallest value given the number of randamization - limit_researched (Yes = 1, No = 0)
limit_reached = if_else(abs((1/p_value) - no_randomizations) <= 1, 1, 0))
We first conducted a series of meta-regression models with each of the moderators introduced above.
# meta-regression: mutiple intercepts
mr_symbiosis1 <- rma.mv(yi = Zr, V = VZr, mods = ~symbiosis - 1, test = "t", random = ~1 |
authors, data = dat)
# meta-regression: contrast
mr_symbiosis2 <- rma.mv(yi = Zr, V = VZr, mods = ~symbiosis, test = "t", random = ~1 |
authors, data = dat)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (Nakagawa & Schielzeth 2013) (R2) from the meta-regression with symbiosis. Note that mu means the group mean while beta represents the contrast between two groups in the Unit column.
# getting marginal R2
r2_symbiosis1 <- R2(mr_symbiosis1)
# getting estimates``
res_symbiosis1 <- get_est(mr_symbiosis1, mod = "symbiosis")
res_symbiosis2 <- get_est(mr_symbiosis2, mod = "symbiosis")
# creating a table
t_symbiosis <- tibble(`Fixed effect` = c(rep(as.character(res_symbiosis1$name), 2),
cont_gen(res_symbiosis1$name)), Unit = c(rep(c("Zr (mu)", "r (mu)"), each = 2),
"Zr (beta)"), Estimate = c(res_symbiosis1$estimate, tanh(res_symbiosis1$estimate),
res_symbiosis2$estimate[2]), `Lower CI [0.025]` = c(res_symbiosis1$lowerCL, tanh(res_symbiosis1$lowerCL),
res_symbiosis2$lowerCL[2]), `Upper CI [0.975]` = c(res_symbiosis1$upperCL, tanh(res_symbiosis1$upperCL),
res_symbiosis2$upperCL[2]), `V[authors]` = c(mr_symbiosis1$sigma2, rep(NA, 4)),
R2 = c(r2_symbiosis1[1], rep(NA, 4)))
t_symbiosis %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| Mutualist | Zr (mu) | 0.631 | 0.534 | 0.727 | 0.078 | 0.054 |
| Parasite | Zr (mu) | 0.487 | 0.419 | 0.554 | NA | NA |
| Mutualist | r (mu) | 0.559 | 0.488 | 0.621 | NA | NA |
| Parasite | r (mu) | 0.451 | 0.396 | 0.503 | NA | NA |
| Mutualist-Parasite | Zr (beta) | -0.144 | -0.260 | -0.028 | NA | NA |
# adding sample size (k) for each category
k_symbiosis <- dat %>% group_by(symbiosis) %>% count()
# getting estimates and predicitons
pred_symbiosis <- get_pred(mr_symbiosis1, mod = "symbiosis")
res_symbiosis1 <- left_join(res_symbiosis1, k_symbiosis, by = c("name" = "symbiosis")) %>% left_join(pred_symbiosis)
#res_symbiosis1
# drawing a funnel plot - fig 2b
fig_symbiosis <- ggplot(data = res_symbiosis1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(symbiosis)),
aes(x= tanh(Zr), y = symbiosis, size = ((1/VZr) + 3), colour = symbiosis), groupOnX = FALSE, alpha=0.2) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("Mutualist" = "#E69F00", "Parasite" = "#56B4E9")) +
scale_fill_manual(values = c("Mutualist" = "#E69F00", "Parasite" = "#56B4E9")) +
annotate('text', x = 0.93, y = c(1.15, 2.15), label= paste("italic(k)==", res_symbiosis1$n), parse = TRUE, hjust = "left", size = 3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")) ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position= c(0, 1), legend.justification = c(0,1)) +
theme(legend.direction="horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_mutualism), xmin = -1, xmax = -0.8, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -1, xmax = -0.8, ymin = 1.6, ymax = 2.2)
fig_symbiosis
# fig 3
b <- ggplot(data = res_symbiosis1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(symbiosis)),
aes(x= tanh(Zr), y = symbiosis, size = ((1/VZr) + 3), colour = symbiosis), groupOnX = FALSE, alpha=0.2) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("Mutualist" = "#E69F00", "Parasite" = "#56B4E9")) +
scale_fill_manual(values = c("Mutualist" = "#E69F00", "Parasite" = "#56B4E9")) +
annotate('text', x = 0.93, y = c(1.15, 2.15), label= paste("italic(k)==", res_symbiosis1$n), parse = TRUE, hjust = "left", size = 3.5) +
labs(x = "", y = "", tag = "b") +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position="none") +
theme(axis.text.y = element_text(size = 10, colour ="black",hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_mutualism), xmin = -1, xmax = -0.8, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -1, xmax = -0.8, ymin = 1.6, ymax = 2.2)
Figure 3b: An orchard plot showing the group-wise means (the categorical variable symbiosis) with their 95% confidence intervals (thick lines) and 95% prediction intervals (thin lines), with observed effect sizes based on various sample sizes.
# reordering
dat$host_tax_broad <- factor(dat$host_tax_broad, levels = c("Microbe", "Plant", "Invert",
"Vert"))
# meta-regression: mutiple intercepts
mr_host_tax_broad1 <- rma.mv(yi = Zr, V = VZr, mods = ~host_tax_broad - 1, test = "t",
random = ~1 | authors, data = dat)
# meta-regression: contrast 1
mr_host_tax_broad2 <- rma.mv(yi = Zr, V = VZr, mods = ~host_tax_broad, test = "t",
random = ~1 | authors, data = dat)
# meta-regression: contrast 2
mr_host_tax_broad3 <- rma.mv(yi = Zr, V = VZr, mods = ~relevel(host_tax_broad, ref = "Plant"),
test = "t", random = ~1 | authors, data = dat)
# meta-regression: contrast 3
mr_host_tax_broad4 <- rma.mv(yi = Zr, V = VZr, mods = ~relevel(host_tax_broad, ref = "Invert"),
test = "t", random = ~1 | authors, data = dat)
Regression coefficients (estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with host_tax_broad. Note that mu means the group mean while beta represents the contrast between two groups in the Unit column.
# getting marginal R2
r2_host_tax_broad1 <- R2(mr_host_tax_broad1)
# getting estimates
res_host_tax_broad1 <- get_est(mr_host_tax_broad1, mod = "host_tax_broad")
res_host_tax_broad2 <- get_est(mr_host_tax_broad2, mod = "host_tax_broad")
# the name bit does not work if relevel....
res_host_tax_broad3 <- get_est(mr_host_tax_broad3, mod = "host_tax_broad")
res_host_tax_broad4 <- get_est(mr_host_tax_broad4, mod = "host_tax_broad")
# creating a table
t_host_tax <- tibble(`Fixed effect` = c(rep(as.character(res_host_tax_broad1$name),
2), cont_gen(res_host_tax_broad1$name)), Unit = c(rep(c("Zr (mu)", "r (mu)"),
each = 4), rep("Zr (beta)", 6)), Estimate = c(res_host_tax_broad1$estimate, tanh(res_host_tax_broad1$estimate),
res_host_tax_broad2$estimate[-1], res_host_tax_broad3$estimate[-(1:2)], res_host_tax_broad4$estimate[-(1:3)]),
`Lower CI [0.025]` = c(res_host_tax_broad1$lowerCL, tanh(res_host_tax_broad1$lowerCL),
res_host_tax_broad2$lowerCL[-1], res_host_tax_broad3$lowerCL[-(1:2)], res_host_tax_broad4$lowerCL[-(1:3)]),
`Upper CI [0.975]` = c(res_host_tax_broad1$upperCL, tanh(res_host_tax_broad1$upperCL),
res_host_tax_broad2$upperCL[-1], res_host_tax_broad3$upperCL[-(1:2)], res_host_tax_broad4$upperCL[-(1:3)]),
`V[authors]` = c(mr_host_tax_broad1$sigma2, rep(NA, 13)), R2 = c(r2_host_tax_broad1[1],
rep(NA, 13)))
t_host_tax %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left") %>%
scroll_box(width = "100%", height = "300px")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| Microbe | Zr (mu) | 0.883 | 0.639 | 1.126 | 0.077 | 0.162 |
| Plant | Zr (mu) | 0.377 | 0.256 | 0.498 | NA | NA |
| Invert | Zr (mu) | 0.629 | 0.522 | 0.736 | NA | NA |
| Vert | Zr (mu) | 0.515 | 0.435 | 0.595 | NA | NA |
| Microbe | r (mu) | 0.708 | 0.564 | 0.810 | NA | NA |
| Plant | r (mu) | 0.360 | 0.250 | 0.460 | NA | NA |
| Invert | r (mu) | 0.557 | 0.479 | 0.627 | NA | NA |
| Vert | r (mu) | 0.474 | 0.410 | 0.534 | NA | NA |
| Microbe-Plant | Zr (beta) | -0.506 | -0.778 | -0.234 | NA | NA |
| Microbe-Invert | Zr (beta) | -0.254 | -0.519 | 0.012 | NA | NA |
| Microbe-Vert | Zr (beta) | -0.367 | -0.624 | -0.111 | NA | NA |
| Plant-Invert | Zr (beta) | 0.252 | 0.091 | 0.414 | NA | NA |
| Plant-Vert | Zr (beta) | 0.139 | -0.006 | 0.284 | NA | NA |
| Invert-Vert | Zr (beta) | -0.114 | -0.244 | 0.017 | NA | NA |
# getting images
image_invertebrate_host <- readPNG(here("images/invertebrate_host_transparentbg.png"))
image_microbe_host <- readPNG(here("images/microbe_host_transparentbg.png"))
image_vertebrate_host <- readPNG(here("images/vertebrate_host_transparentbg.png"))
image_plant_host <- readPNG(here("images/plant_host_transparentbg.png"))
# adding sample size (k) for each category
k_host_tax_broad <- dat %>% group_by(host_tax_broad) %>% count()
# getting estimates and predicitons
pred_host_tax_broad <- get_pred(mr_host_tax_broad1, mod = "host_tax_broad")
res_host_tax_broad1 <- left_join(res_host_tax_broad1, k_host_tax_broad, by = c("name" = "host_tax_broad")) %>% left_join(pred_host_tax_broad)
#res_symbiosis1
# drawing a funnel plot - fig 2b
fig_host_tax_broad <- ggplot(data = res_host_tax_broad1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(host_tax_broad)),
aes(x= tanh(Zr), y = host_tax_broad, size = ((1/VZr) + 3), colour = host_tax_broad), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("Microbe" = "#009E73", "Plant" = "#F0E422", "Invert"= "#0072B2", "Vert" = "#D55E00")) +
scale_fill_manual(values = c("Microbe" = "#009E73", "Plant" = "#F0E422", "Invert"= "#0072B2", "Vert" = "#D55E00")) +
scale_y_discrete(labels = c("Microbe" = "Microbe", "Plant" = "Plant", "Invert"= "Invertebrate", "Vert" = "Vertebrate")) +
annotate('text', x = 0.93, y = 1:4 + 0.15, label= paste("italic(k)==", res_host_tax_broad1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")) ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position= c(0, 1), legend.justification = c(0,1)) +
theme(legend.direction="horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_microbe_host), xmin = -1, xmax = -0.8, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_plant_host), xmin = -1, xmax = -0.8, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_invertebrate_host), xmin = -1, xmax = -0.8, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_vertebrate_host), xmin = -1, xmax = -0.8, ymin = 3.6, ymax = 4.2)
fig_host_tax_broad
# fig 3c
c <- ggplot(data = res_host_tax_broad1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(host_tax_broad)),
aes(x= tanh(Zr), y = host_tax_broad, size = ((1/VZr) + 3), colour = host_tax_broad), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("Microbe" = "#009E73", "Plant" = "#F0E422", "Invert"= "#0072B2", "Vert" = "#D55E00")) +
scale_fill_manual(values = c("Microbe" = "#009E73", "Plant" = "#F0E422", "Invert"= "#0072B2", "Vert" = "#D55E00")) +
scale_y_discrete(labels = c("Microbe" = "Microbe", "Plant" = "Plant", "Invert"= "Invertebrate", "Vert" = "Vertebrate")) +
annotate('text', x = 0.93, y = 1:4 + 0.15, label= paste("italic(k)==", res_host_tax_broad1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = "", y = "", size = expression(paste(italic(n), " (# of species pairs)")) , tag = "c") +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position="none") +
theme(axis.text.y = element_text(size = 10, colour ="black",hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_microbe_host), xmin = -1, xmax = -0.8, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_plant_host), xmin = -1, xmax = -0.8, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_invertebrate_host), xmin = -1, xmax = -0.8, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_vertebrate_host), xmin = -1, xmax = -0.8, ymin = 3.6, ymax = 4.2)
Figure 3c: An orchard plot showing the group-wise means (the categorical variable host_tax_broad) with their 95% confidence intervals (thick lines) and 95% prediction intervals (thin lines), with observed effect sizes based on various sample sizes.
# reordering
dat$symbiont_tax_broad <- factor(dat$symbiont_tax_broad, levels = c("Microbe", "Plant",
"Invert", "Vert"))
# sizes <- factor(sizes, levels = c('small', 'medium', 'large')) sizes > [1]
# small large large small medium > Levels: small medium large meta-regression:
# mutiple intercepts
mr_symbiont_tax_broad1 <- rma.mv(yi = Zr, V = VZr, mods = ~symbiont_tax_broad - 1,
test = "t", random = ~1 | authors, data = dat)
# meta-regression: contrast 1
mr_symbiont_tax_broad2 <- rma.mv(yi = Zr, V = VZr, mods = ~symbiont_tax_broad, test = "t",
random = ~1 | authors, data = dat)
# meta-regression: contrast 2
mr_symbiont_tax_broad3 <- rma.mv(yi = Zr, V = VZr, mods = ~relevel(symbiont_tax_broad,
ref = "Plant"), test = "t", random = ~1 | authors, data = dat)
# meta-regression: contrast 3
mr_symbiont_tax_broad4 <- rma.mv(yi = Zr, V = VZr, mods = ~relevel(symbiont_tax_broad,
ref = "Invert"), test = "t", random = ~1 | authors, data = dat)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with symbiont_tax_broad. Note that mu means the group mean while beta represents the contrast between two groups in the Unit column.
# getting marginal R2
r2_symbiont_tax_broad1 <- R2(mr_symbiont_tax_broad1)
# getting estimates
res_symbiont_tax_broad1 <- get_est(mr_symbiont_tax_broad1, mod = "symbiont_tax_broad")
res_symbiont_tax_broad2 <- get_est(mr_symbiont_tax_broad2, mod = "symbiont_tax_broad")
res_symbiont_tax_broad3 <- get_est(mr_symbiont_tax_broad3, mod = "symbiont_tax_broad")
res_symbiont_tax_broad4 <- get_est(mr_symbiont_tax_broad4, mod = "symbiont_tax_broad")
# creating a table
tibble(`Fixed effect` = c(as.character(res_symbiont_tax_broad1$name), as.character(res_symbiont_tax_broad1$name),
cont_gen(res_symbiont_tax_broad1$name)), Unit = c(rep(c("Zr (mu)", "r (mu)"),
each = 4), rep("Zr (beta)", 6)), Estimate = c(res_symbiont_tax_broad1$estimate,
tanh(res_symbiont_tax_broad1$estimate), res_symbiont_tax_broad2$estimate[-1],
res_symbiont_tax_broad3$estimate[-(1:2)], res_symbiont_tax_broad4$estimate[-(1:3)]),
`Lower CI [0.025]` = c(res_symbiont_tax_broad1$lowerCL, tanh(res_symbiont_tax_broad1$lowerCL),
res_symbiont_tax_broad2$lowerCL[-1], res_symbiont_tax_broad3$lowerCL[-(1:2)],
res_symbiont_tax_broad4$lowerCL[-(1:3)]), `Upper CI [0.975]` = c(res_symbiont_tax_broad1$upperCL,
tanh(res_symbiont_tax_broad1$upperCL), res_symbiont_tax_broad2$upperCL[-1],
res_symbiont_tax_broad3$upperCL[-(1:2)], res_symbiont_tax_broad4$upperCL[-(1:3)]),
`V[authors]` = c(mr_symbiont_tax_broad1$sigma2, rep(NA, 13)), R2 = c(r2_symbiont_tax_broad1[1],
rep(NA, 13))) %>% kable("html", digits = 3) %>% kable_styling("striped",
position = "left") %>% scroll_box(width = "100%", height = "300px")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| Microbe | Zr (mu) | 0.546 | 0.467 | 0.624 | 0.081 | 0.078 |
| Plant | Zr (mu) | 1.140 | 0.671 | 1.609 | NA | NA |
| Invert | Zr (mu) | 0.507 | 0.423 | 0.591 | NA | NA |
| Vert | Zr (mu) | 0.440 | -0.191 | 1.071 | NA | NA |
| Microbe | r (mu) | 0.497 | 0.436 | 0.554 | NA | NA |
| Plant | r (mu) | 0.814 | 0.585 | 0.923 | NA | NA |
| Invert | r (mu) | 0.467 | 0.399 | 0.531 | NA | NA |
| Vert | r (mu) | 0.414 | -0.189 | 0.790 | NA | NA |
| Microbe-Plant | Zr (beta) | 0.594 | 0.118 | 1.070 | NA | NA |
| Microbe-Invert | Zr (beta) | -0.039 | -0.154 | 0.076 | NA | NA |
| Microbe-Vert | Zr (beta) | -0.106 | -0.742 | 0.531 | NA | NA |
| Plant-Invert | Zr (beta) | -0.633 | -1.110 | -0.156 | NA | NA |
| Plant-Vert | Zr (beta) | -0.700 | -1.487 | 0.087 | NA | NA |
| Invert-Vert | Zr (beta) | -0.067 | -0.704 | 0.570 | NA | NA |
# getting images
image_invertebrate_parasite <- readPNG(here("images/invertebrate_parasite_transparentbg.png"))
image_microbe_parasite <- readPNG(here("images/microbe_parasite_transparentbg.png"))
image_vertebrate_parasite <- readPNG(here("images/vertebrate_parasite_transparentbg.png"))
image_plant_parasite <- readPNG(here("images/plant_parasite_transparentbg.png"))
# adding sample size (k) for each category
k_symbiont_tax_broad <- dat %>% group_by(symbiont_tax_broad) %>% count()
# getting estimates and predicitons
pred_symbiont_tax_broad <- get_pred(mr_symbiont_tax_broad1, mod = "symbiont_tax_broad")
res_symbiont_tax_broad1 <- left_join(res_symbiont_tax_broad1, k_symbiont_tax_broad, by = c("name" = "symbiont_tax_broad")) %>% left_join(pred_symbiont_tax_broad)
#res_symbiosis1
# drawing a funnel plot - fig 2b
fig_symbiont_tax_broad <- ggplot(data = res_symbiont_tax_broad1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(symbiont_tax_broad)),
aes(x= tanh(Zr), y = symbiont_tax_broad, size = ((1/VZr) + 3), colour = symbiont_tax_broad), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("Microbe" = "#009E73", "Plant" = "#F0E422", "Invert"= "#0072B2", "Vert" = "#D55E00" )) +
scale_fill_manual(values = c("Microbe" = "#009E73", "Plant" = "#F0E422", "Invert"= "#0072B2", "Vert" = "#D55E00")) +
scale_y_discrete(labels = c("Microbe" = "Microbe", "Plant" = "Plant", "Invert"= "Invertebrate", "Vert" = "Vertebrate")) +
annotate('text', x = 0.93, y = 1:4 + 0.15, label= paste("italic(k)==", res_symbiont_tax_broad1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")) ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position= c(0, 1), legend.justification = c(0,1)) +
theme(legend.direction="horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_microbe_parasite), xmin = -1, xmax = -0.8, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_plant_parasite), xmin = -1, xmax = -0.8, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_invertebrate_parasite), xmin = -1, xmax = -0.8, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_vertebrate_parasite), xmin = -1, xmax = -0.8, ymin = 3.6, ymax = 4.2)
fig_symbiont_tax_broad
# fig 3d
d <- ggplot(data = res_symbiont_tax_broad1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(symbiont_tax_broad)),
aes(x= tanh(Zr), y = symbiont_tax_broad, size = ((1/VZr) + 3), colour = symbiont_tax_broad), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("Microbe" = "#009E73", "Plant" = "#F0E422", "Invert"= "#0072B2", "Vert" = "#D55E00" )) +
scale_fill_manual(values = c("Microbe" = "#009E73", "Plant" = "#F0E422", "Invert"= "#0072B2", "Vert" = "#D55E00")) +
scale_y_discrete(labels = c("Microbe" = "Microbe", "Plant" = "Plant", "Invert"= "Invertebrate", "Vert" = "Vertebrate")) +
annotate('text', x = 0.93, y = 1:4 + 0.15, label= paste("italic(k)==", res_symbiont_tax_broad1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")), tag = "d") +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position="none") +
theme(axis.text.y = element_text(size = 10, colour ="black",hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_microbe_parasite), xmin = -1, xmax = -0.8, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_plant_parasite), xmin = -1, xmax = -0.8, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_invertebrate_parasite), xmin = -1, xmax = -0.8, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_vertebrate_parasite), xmin = -1, xmax = -0.8, ymin = 3.6, ymax = 4.2)
Figure 3d: An orchard plot showing the group-wise means (the categorical variable symbiont_tax_broad) with their 95% confidence intervals (thick lines) and 95% prediction intervals (thin lines), with observed effect sizes based on various sample sizes.
# meta-regression
mr_host_range_link_ratio <- rma.mv(yi = Zr, V = VZr, mods = ~log(host_range_link_ratio),
random = ~1 | authors, data = dat)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with log(host_range_link_ratio). Note that mu shows the overall while beta represents a slope in the Unit column.
# getting marginal R2
r2_host_range_link_ratio <- R2(mr_host_range_link_ratio)
# getting estimates: name does not work for slopes
res_host_range_link_ratio <- get_est(mr_host_range_link_ratio, mod = "log(host_range_link_ratio)")
# creating a table
tibble(`Fixed effect` = c("Intercept", "Intercept", "log(host_range_link_ratio)"),
Unit = c(rep(c("Zr (mu)", "r (mu)"), 1), rep("Zr (beta)", 1)), Estimate = c(res_host_range_link_ratio$estimate[1],
tanh(res_host_range_link_ratio$estimate[1]), res_host_range_link_ratio$estimate[2]),
`Lower CI [0.025]` = c(res_host_range_link_ratio$lowerCL[1], tanh(res_host_range_link_ratio$lowerCL[1]),
res_host_range_link_ratio$lowerCL[2]), `Upper CI [0.975]` = c(res_host_range_link_ratio$upperCL[1],
tanh(res_host_range_link_ratio$upperCL[1]), res_host_range_link_ratio$upperCL[2]),
`V[authors]` = c(mr_host_range_link_ratio$sigma2, NA, NA), R2 = c(r2_host_range_link_ratio[1],
NA, NA)) %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| Intercept | Zr (mu) | 0.549 | 0.484 | 0.614 | 0.083 | 0.002 |
| Intercept | r (mu) | 0.500 | 0.449 | 0.547 | NA | NA |
| log(host_range_link_ratio) | Zr (beta) | -0.037 | -0.187 | 0.114 | NA | NA |
# newmods <- seq(-0.3, 2.2, by = 0.1)
# pred_host_range_link_ratio <-predict.rma(mr_host_range_link_ratio, newmods = newmods)
# ribbon_dat <- tibble(newmods = newmods, ymin = pred_host_range_link_ratio$ci.lb, ymax = pred_host_range_link_ratio$ci.ub)
pred_host_range_link_ratio <-predict.rma(mr_host_range_link_ratio)
# plotting
fig_host_range_link_ratio <- dat %>%
filter(!is.na(host_range_link_ratio)) %>% # getting ride of NA values
mutate(ymin = pred_host_range_link_ratio$ci.lb,
ymax = pred_host_range_link_ratio$ci.ub,
ymin2 = pred_host_range_link_ratio$cr.lb,
ymax2 = pred_host_range_link_ratio$cr.ub,
pred = pred_host_range_link_ratio$pred) %>%
ggplot(aes(x = log(host_range_link_ratio), y = Zr, size = (1/VZr) + 3, )) +
geom_point(shape = 21, fill = "grey90") +
#geom_ribbon(aes(ymin = ymin, ymax = ymax), fill = "#0072B2") + # not quite sure why this does not work
geom_smooth(aes(y = ymin2), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25, colour = "#0072B2") +
geom_smooth(aes(y = ymax2), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25, colour = "#0072B2") +
geom_smooth(aes(y = ymin), method = "loess", se = FALSE,lty = "dotted", lwd = 0.25, colour ="#D55E00") +
geom_smooth(aes(y = ymax), method = "loess", se = FALSE, lty ="dotted", lwd = 0.25, colour ="#D55E00") +
geom_smooth(aes(y = pred), method = "loess", se = FALSE, lty ="dashed", lwd = 0.5, colour ="black") +
ylim(-1, 2) + xlim(-0.05, 2) +
#geom_abline(intercept = mr_host_range_link_ratio$beta[[1]], slope = mr_host_range_link_ratio$beta[[2]], alpha = 0.7, linetype = "dashed", size = 0.5) +
labs(x = "ln(host range link ratio)", y = expression(paste(italic(Zr), " (effect size)")), size = expression(paste(italic(n), " (# of species pairs)"))) +
guides(fill = "none", colour = "none") +
# themses
theme_bw() +
theme(legend.position= c(1, 1), legend.justification = c(1, 1)) +
theme(legend.direction="horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90))
fig_host_range_link_ratio
A bubble plot showing a predicted regression line for the continuous variable log(host_range_link_ratio), indicating 95% confidence regions (orange dotted lines) and 95% prediction regions (blue dotted lines) with observed effect sizes based on various sample sizes.
# meta-regression
mr_host_range_taxonomic_breadth <- rma.mv(yi = Zr, V = VZr, mods = ~log(host_range_taxonomic_breadth),
random = ~1 | authors, data = dat)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with log(host_range_taxonomic_breadth). Note that mu shows the overall while beta represents a slope in the Unit column.
# getting marginal R2
r2_host_range_taxonomic_breadth <- R2(mr_host_range_taxonomic_breadth)
# getting estimates: name does not work for slopes
res_host_range_taxonomic_breadth <- get_est(mr_host_range_taxonomic_breadth, mod = "log(host_range_taxonomic_breadth)")
# creating a table
tibble(`Fixed effect` = c("Intercept", "Intercept", "log(host_range_taxonomic_breadth)"),
Unit = c(rep(c("Zr (mu)", "r (mu)"), 1), rep("Zr (beta)", 1)), Estimate = c(res_host_range_taxonomic_breadth$estimate[1],
tanh(res_host_range_taxonomic_breadth$estimate[1]), res_host_range_taxonomic_breadth$estimate[2]),
`Lower CI [0.025]` = c(res_host_range_taxonomic_breadth$lowerCL[1], tanh(res_host_range_taxonomic_breadth$lowerCL[1]),
res_host_range_taxonomic_breadth$lowerCL[2]), `Upper CI [0.975]` = c(res_host_range_taxonomic_breadth$upperCL[1],
tanh(res_host_range_taxonomic_breadth$upperCL[1]), res_host_range_taxonomic_breadth$upperCL[2]),
`V[authors]` = c(mr_host_range_taxonomic_breadth$sigma2, NA, NA), R2 = c(r2_host_range_taxonomic_breadth[1],
NA, NA)) %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| Intercept | Zr (mu) | 0.559 | 0.494 | 0.623 | 0.081 | 0.001 |
| Intercept | r (mu) | 0.507 | 0.458 | 0.554 | NA | NA |
| log(host_range_taxonomic_breadth) | Zr (beta) | -0.036 | -0.216 | 0.143 | NA | NA |
pred_host_range_taxonomic_breadth <-predict.rma(mr_host_range_taxonomic_breadth)
# plotting
fig_host_range_taxonomic_breadth <- dat %>%
filter(!is.na(host_range_taxonomic_breadth)) %>% # getting ride of NA values
mutate(ymin = pred_host_range_taxonomic_breadth$ci.lb,
ymax = pred_host_range_taxonomic_breadth$ci.ub,
ymin2 = pred_host_range_taxonomic_breadth$cr.lb,
ymax2 = pred_host_range_taxonomic_breadth$cr.ub,
pred = pred_host_range_taxonomic_breadth$pred) %>%
ggplot(aes(x = log(host_range_taxonomic_breadth), y = Zr, size = (1/VZr) + 3, )) +
geom_point(shape = 21, fill = "grey90") +
#geom_ribbon(aes(ymin = ymin, ymax = ymax), fill = "#0072B2") + # not quite sure why this does not work
geom_smooth(aes(y = ymin2), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25, colour = "#0072B2") +
geom_smooth(aes(y = ymax2), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25, colour = "#0072B2") +
geom_smooth(aes(y = ymin), method = "loess", se = FALSE,lty = "dotted", lwd = 0.25, colour ="#D55E00") +
geom_smooth(aes(y = ymax), method = "loess", se = FALSE, lty ="dotted", lwd = 0.25, colour ="#D55E00") +
geom_smooth(aes(y = pred), method = "loess", se = FALSE, lty ="dashed", lwd = 0.5, colour ="black") +
ylim(-1, 2) + xlim(0, 1.5) +
#geom_abline(intercept = mr_host_range_link_ratio$beta[[1]], slope = mr_host_range_link_ratio$beta[[2]], alpha = 0.7, linetype = "dashed", size = 0.5) +
labs(x = "ln(host range taxonomic breadth)", y = expression(paste(italic(Zr), " (effect size)")), size = expression(paste(italic(n), " (# of species pairs)"))) +
guides(fill = "none", colour = "none") +
# themses
theme_bw() +
theme(legend.position= c(1, 1), legend.justification = c(1, 1)) +
theme(legend.direction="horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90))
fig_host_range_taxonomic_breadth
A bubble plot showing a predicted regression line for the continuous variable log(log(host_range_taxonomic_breadth), indicating 95% confidence regions (orange dotted lines) and 95% prediction regions (blue dotted lines) with observed effect sizes based on various sample sizes.
# reordering
dat$endo_or_ecto <- factor(dat$endo_or_ecto, levels = c("Endo/Ecto", "Endo", "Ecto"))
# meta-regression: mutiple intercepts
mr_endo_or_ecto1 <- rma.mv(yi = Zr, V = VZr, mods = ~endo_or_ecto - 1, test = "t",
random = ~1 | authors, data = dat)
# meta-regression: contrast 1
mr_endo_or_ecto2 <- rma.mv(yi = Zr, V = VZr, mods = ~endo_or_ecto, test = "t", random = ~1 |
authors, data = dat)
# meta-regression: contrast 2
mr_endo_or_ecto3 <- rma.mv(yi = Zr, V = VZr, mods = ~relevel(endo_or_ecto, ref = "Endo"),
test = "t", random = ~1 | authors, data = dat)
Regression coefficients (estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with endo_or_ecto. Note that mu means the group mean while beta represents the contrast between two groups in the Unit column.
# getting marginal R2
r2_endo_or_ecto1 <- R2(mr_endo_or_ecto1)
# getting estimates
res_endo_or_ecto1 <- get_est(mr_endo_or_ecto1, mod = "endo_or_ecto")
res_endo_or_ecto2 <- get_est(mr_endo_or_ecto2, mod = "endo_or_ecto")
res_endo_or_ecto3 <- get_est(mr_endo_or_ecto3, mod = "endo_or_ecto")
# creating a table
tibble(`Fixed effect` = c(as.character(res_endo_or_ecto1$name), as.character(res_endo_or_ecto1$name),
cont_gen(res_endo_or_ecto1$name)), Unit = c(rep(c("Zr (mu)", "r (mu)"), each = 3),
rep("Zr (beta)", 3)), Estimate = c(res_endo_or_ecto1$estimate, tanh(res_endo_or_ecto1$estimate),
res_endo_or_ecto2$estimate[-1], res_endo_or_ecto3$estimate[-(1:2)]), `Lower CI [0.025]` = c(res_endo_or_ecto1$lowerCL,
tanh(res_endo_or_ecto1$lowerCL), res_endo_or_ecto2$lowerCL[-1], res_endo_or_ecto3$lowerCL[-(1:2)]),
`Upper CI [0.975]` = c(res_endo_or_ecto1$upperCL, tanh(res_endo_or_ecto1$upperCL),
res_endo_or_ecto2$upperCL[-1], res_endo_or_ecto3$upperCL[-(1:2)]), `V[authors]` = c(mr_endo_or_ecto1$sigma2,
rep(NA, 8)), R2 = c(r2_endo_or_ecto1[1], rep(NA, 8))) %>% kable("html", digits = 3) %>%
kable_styling("striped", position = "left")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| Endo/Ecto | Zr (mu) | 0.492 | 0.029 | 0.955 | 0.083 | 0.002 |
| Endo | Zr (mu) | 0.530 | 0.460 | 0.599 | NA | NA |
| Ecto | Zr (mu) | 0.555 | 0.451 | 0.659 | NA | NA |
| Endo/Ecto | r (mu) | 0.456 | 0.029 | 0.742 | NA | NA |
| Endo | r (mu) | 0.485 | 0.430 | 0.536 | NA | NA |
| Ecto | r (mu) | 0.504 | 0.423 | 0.578 | NA | NA |
| Endo/Ecto-Endo | Zr (beta) | 0.037 | -0.431 | 0.506 | NA | NA |
| Endo/Ecto-Ecto | Zr (beta) | 0.063 | -0.412 | 0.537 | NA | NA |
| Endo-Ecto | Zr (beta) | 0.026 | -0.100 | 0.151 | NA | NA |
# getting images
image_endoparasite <- readPNG(here("images/endoparasite_transparentbg.png"))
image_ectoparasite <- readPNG(here("images/ectoparasite_transparentbg.png"))
# adding sample size (k) for each category
k_endo_or_ecto <- dat %>% group_by(endo_or_ecto) %>% count()
# getting estimates and predicitons
pred_endo_or_ecto <- get_pred(mr_endo_or_ecto1, mod = "endo_or_ecto")
res_endo_or_ecto1 <- left_join(res_endo_or_ecto1, k_endo_or_ecto, by = c("name" = "endo_or_ecto")) %>% left_join(pred_endo_or_ecto)
#res_symbiosis1
# drawing a funnel plot - fig 2b
fig_endo_or_ecto <- ggplot(data = res_endo_or_ecto1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(endo_or_ecto)),
aes(x= tanh(Zr), y = endo_or_ecto, size = ((1/VZr) + 3), colour = endo_or_ecto), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("Endo/Ecto" = "#0072B2", "Endo" = "#D55E00", "Ecto"= "#CC79A7")) +
scale_fill_manual(values = c("Endo/Ecto" = "#0072B2", "Endo" = "#D55E00", "Ecto"= "#CC79A7")) +
scale_y_discrete(labels = c("Endo/Ecto" = "Both", "Endo" = "Endosymbiosis", "Ecto"= "Ectosymbiosis")) +
annotate('text', x = 0.93, y = 1:3 + 0.15, label= paste("italic(k)==", res_endo_or_ecto1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")) ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position= c(0, 1), legend.justification = c(0,1)) +
theme(legend.direction = "horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90)) +
# adding images
annotation_custom(rasterGrob(image_endoparasite), xmin = -1, xmax = -0.8, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_ectoparasite), xmin = -1, xmax = -0.8, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_ectoparasite), xmin = -1.1, xmax = -0.9, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_endoparasite), xmin = -0.9, xmax = -0.7, ymin = 0.6, ymax = 1.2)
fig_endo_or_ecto
# fig 3e
e <- ggplot(data = res_endo_or_ecto1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(endo_or_ecto)),
aes(x= tanh(Zr), y = endo_or_ecto, size = ((1/VZr) + 3), colour = endo_or_ecto), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("Endo/Ecto" = "#0072B2", "Endo" = "#D55E00", "Ecto"= "#CC79A7")) +
scale_fill_manual(values = c("Endo/Ecto" = "#0072B2", "Endo" = "#D55E00", "Ecto"= "#CC79A7")) +
scale_y_discrete(labels = c("Endo/Ecto" = "Both", "Endo" = "Endosymbiosis", "Ecto"= "Ectosymbiosis")) +
annotate('text', x = 0.93, y = 1:3 + 0.15, label= paste("italic(k)==", res_endo_or_ecto1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = "", y = "", size = expression(paste(italic(n), " (# of species pairs)")), tag = "e" ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position="none") +
theme(axis.text.y = element_text(size = 10, colour ="black",hjust = 0.5, angle = 90)) +
# adding images
annotation_custom(rasterGrob(image_endoparasite), xmin = -1, xmax = -0.8, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_ectoparasite), xmin = -1, xmax = -0.8, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_ectoparasite), xmin = -1.1, xmax = -0.9, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_endoparasite), xmin = -0.9, xmax = -0.7, ymin = 0.6, ymax = 1.2)
Figure 4a: An orchard plot showing the group-wise means (the categorical variable endo_or_ecto) with their 95% confidence intervals (thick lines) and 95% prediction intervals (thin lines), with observed effect sizes based on various sample sizes.
# meta-regression: mutiple intercepts
mr_mode_of_transmission_broad1 <- rma.mv(yi = Zr, V = VZr, mods = ~mode_of_transmission_broad -
1, test = "t", random = ~1 | authors, data = dat)
# meta-regression: contrast 1
mr_mode_of_transmission_broad2 <- rma.mv(yi = Zr, V = VZr, mods = ~mode_of_transmission_broad,
test = "t", random = ~1 | authors, data = dat)
# meta-regression: contrast 2
mr_mode_of_transmission_broad3 <- rma.mv(yi = Zr, V = VZr, mods = ~relevel(mode_of_transmission_broad,
ref = "vertical"), test = "t", random = ~1 | authors, data = dat)
Regression coefficients (estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with mode_of_transmission_broad. Note that mu means the group mean while beta represents the contrast between two groups in the Unit column.
# getting marginal R2
r2_mode_of_transmission_broad1 <- R2(mr_mode_of_transmission_broad1)
# getting estimates
res_mode_of_transmission_broad1 <- get_est(mr_mode_of_transmission_broad1, mod = "mode_of_transmission_broad")
res_mode_of_transmission_broad2 <- get_est(mr_mode_of_transmission_broad2, mod = "mode_of_transmission_broad")
res_mode_of_transmission_broad3 <- get_est(mr_mode_of_transmission_broad3, mod = "mode_of_transmission_broad")
# creating a table
t_transmission <- tibble(`Fixed effect` = c(as.character(res_mode_of_transmission_broad1$name),
as.character(res_mode_of_transmission_broad1$name), cont_gen(res_mode_of_transmission_broad1$name)),
Unit = c(rep(c("Zr (mu)", "r (mu)"), each = 3), rep("Zr (beta)", 3)), Estimate = c(res_mode_of_transmission_broad1$estimate,
tanh(res_mode_of_transmission_broad1$estimate), res_mode_of_transmission_broad2$estimate[-1],
res_mode_of_transmission_broad3$estimate[-(1:2)]), `Lower CI [0.025]` = c(res_mode_of_transmission_broad1$lowerCL,
tanh(res_mode_of_transmission_broad1$lowerCL), res_mode_of_transmission_broad2$lowerCL[-1],
res_mode_of_transmission_broad3$lowerCL[-(1:2)]), `Upper CI [0.975]` = c(res_mode_of_transmission_broad1$upperCL,
tanh(res_mode_of_transmission_broad1$upperCL), res_mode_of_transmission_broad2$upperCL[-1],
res_mode_of_transmission_broad3$upperCL[-(1:2)]), `V[authors]` = c(mr_mode_of_transmission_broad1$sigma2,
rep(NA, 8)), R2 = c(r2_mode_of_transmission_broad1[1], rep(NA, 8)))
t_transmission %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| both | Zr (mu) | 0.578 | 0.463 | 0.693 | 0.061 | 0.187 |
| horizontal | Zr (mu) | 0.446 | 0.378 | 0.515 | NA | NA |
| vertical | Zr (mu) | 0.751 | 0.634 | 0.868 | NA | NA |
| both | r (mu) | 0.521 | 0.432 | 0.600 | NA | NA |
| horizontal | r (mu) | 0.419 | 0.361 | 0.474 | NA | NA |
| vertical | r (mu) | 0.636 | 0.561 | 0.701 | NA | NA |
| both-horizontal | Zr (beta) | -0.131 | -0.265 | 0.002 | NA | NA |
| both-vertical | Zr (beta) | 0.174 | 0.009 | 0.338 | NA | NA |
| horizontal-vertical | Zr (beta) | -0.305 | -0.441 | -0.169 | NA | NA |
# getting images
image_horizontal <- readPNG(here("images/horizontal_transparentbg.png"))
image_vertical <- readPNG(here("images/vertical_transparentbg.png"))
image_both <- readPNG(here("images/horizontal_vertical_transparentbg.png"))
# adding sample size (k) for each category
k_mode_of_transmission_broad <- dat %>% group_by(mode_of_transmission_broad) %>% count()
# getting estimates and predicitons
pred_mode_of_transmission_broad <- get_pred(mr_mode_of_transmission_broad1, mod = "mode_of_transmission_broad")
res_mode_of_transmission_broad1 <- left_join(res_mode_of_transmission_broad1, k_mode_of_transmission_broad, by = c("name" = "mode_of_transmission_broad")) %>% left_join(pred_mode_of_transmission_broad)
#res_symbiosis1
# drawing a funnel plot - fig 2b
fig_mode_of_transmission_broad <- ggplot(data = res_mode_of_transmission_broad1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(mode_of_transmission_broad)),
aes(x= tanh(Zr), y = mode_of_transmission_broad, size = ((1/VZr) + 3), colour = mode_of_transmission_broad), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("both" = "#0072B2", "horizontal" = "#D55E00", "vertical"= "#CC79A7")) +
scale_fill_manual(values = c("both" = "#0072B2", "horizontal" = "#D55E00", "vertical"= "#CC79A7")) +
scale_y_discrete(labels = c("both" = "Both", "horizontal" = "Horizontal", "vertical"= "Vertical")) +
annotate('text', x = 0.93, y = (1:3 + 0.15), label= paste("italic(k)==", res_mode_of_transmission_broad1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")) ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position= c(0, 1), legend.justification = c(0,1)) +
theme(legend.direction = "horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90)) +
# adding images
annotation_custom(rasterGrob(image_horizontal), xmin = -1, xmax = -0.7, ymin = 1.4, ymax = 2.2) +
annotation_custom(rasterGrob(image_vertical), xmin = -1, xmax = -0.7, ymin = 2.4, ymax = 3.2) +
annotation_custom(rasterGrob(image_both), xmin = -1, xmax = -0.7, ymin = 0.4, ymax = 1.2)
fig_mode_of_transmission_broad
# fig 3f
f <- ggplot(data = res_mode_of_transmission_broad1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(mode_of_transmission_broad)),
aes(x= tanh(Zr), y = mode_of_transmission_broad, size = ((1/VZr) + 3), colour = mode_of_transmission_broad), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("both" = "#0072B2", "horizontal" = "#D55E00", "vertical"= "#CC79A7")) +
scale_fill_manual(values = c("both" = "#0072B2", "horizontal" = "#D55E00", "vertical"= "#CC79A7")) +
scale_y_discrete(labels = c("both" = "Both", "horizontal" = "Horizontal", "vertical"= "Vertical")) +
annotate('text', x = 0.93, y = (1:3 + 0.15), label= paste("italic(k)==", res_mode_of_transmission_broad1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = "", y = "", size = expression(paste(italic(n), " (# of species pairs)")), tag = "f" ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position="none") +
theme(axis.text.y = element_text(size = 10, colour ="black",hjust = 0.5, angle = 90)) +
# adding images
annotation_custom(rasterGrob(image_horizontal), xmin = -1, xmax = -0.7, ymin = 1.4, ymax = 2.2) +
annotation_custom(rasterGrob(image_vertical), xmin = -1, xmax = -0.7, ymin = 2.4, ymax = 3.2) +
annotation_custom(rasterGrob(image_both), xmin = -1, xmax = -0.7, ymin = 0.4, ymax = 1.2)
Figure 4b: An orchard plot showing the group-wise means (the categorical variable mode_of_transmission_broad), indicating 95% confidence intervals (thick lines) and 95% prediction intervals (thin lines), with observed effect sizes based on various sample sizes.
# reordering
dat$symbiosis_transmission <- factor(dat$symbiosis_transmission, levels = c("Mutualistboth",
"Mutualisthorizontal", "Mutualistvertical", "Parasiteboth", "Parasitehorizontal"),
labels = c("MutualistBoth", "MutualistHorizontal", "MutualistVertical", "ParasiteBoth",
"ParasiteHorizontal"))
# meta-regression: mutiple intercepts
mr_symbiosis_transmission1 <- rma.mv(yi = Zr, V = VZr, mods = ~symbiosis_transmission -
1, test = "t", random = ~1 | authors, data = dat)
# # meta-regression: contrasts x 10 getting the level names out
level_names <- levels(dat$symbiosis_transmission)
# helper function to run metafor meta-regression
run_rma <- function(name) {
rma.mv(yi = Zr, V = VZr, mods = ~relevel(symbiosis_transmission, ref = name),
test = "t", random = ~1 | authors, data = dat)
}
# results of meta-regression including all contrast results; taking the last
# level out ([-length(level_names)])
mr_symbiosis_transmission <- map(level_names[-length(level_names)], run_rma)
Regression coefficients (estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with symbiosis_transmission. Note that mu means the group mean while beta represents the contrast between two groups in the Unit column.
# getting marginal R2
r2_symbiosis_transmission1 <- R2(mr_symbiosis_transmission1)
# getting estimates
res_symbiosis_transmission1 <- get_est(mr_symbiosis_transmission1, mod = "symbiosis_transmission")
res_symbiosis_transmission <- map(mr_symbiosis_transmission, ~get_est(.x, mod = "symbiosis_transmission"))
# a list of the numbers to take out unnecessary contrasts
contra_list <- Map(seq, from = 1, to = 1:4)
# you need to flatten twice: first to make it a list and make it a vector
estimates <- map2(res_symbiosis_transmission, contra_list, ~.x[-(.y), "estimate"]) %>%
flatten() %>% flatten_dbl()
lowerCLs <- map2(res_symbiosis_transmission, contra_list, ~.x[-(.y), "lowerCL"]) %>%
flatten() %>% flatten_dbl()
upperCLs <- map2(res_symbiosis_transmission, contra_list, ~.x[-(.y), "upperCL"]) %>%
flatten() %>% flatten_dbl()
# creating a table
tibble(`Fixed effect` = c(as.character(res_symbiosis_transmission1$name), as.character(res_symbiosis_transmission1$name),
cont_gen(res_symbiosis_transmission1$name)), Unit = c(rep(c("Zr (mu)", "r (mu)"),
each = 5), rep("Zr (beta)", 10)), Estimate = c(res_symbiosis_transmission1$estimate,
tanh(res_symbiosis_transmission1$estimate), estimates), `Lower CI [0.025]` = c(res_symbiosis_transmission1$lowerCL,
tanh(res_symbiosis_transmission1$lowerCL), lowerCLs), `Upper CI [0.975]` = c(res_symbiosis_transmission1$upperCL,
tanh(res_symbiosis_transmission1$upperCL), upperCLs), `V[authors]` = c(mr_symbiosis_transmission1$sigma2,
rep(NA, (5 + 5 + choose(5, 2)) - 1)), R2 = c(r2_symbiosis_transmission1[1], rep(NA,
(5 + 5 + choose(5, 2)) - 1))) %>% kable("html", digits = 3) %>% kable_styling("striped",
position = "left") %>% scroll_box(width = "100%", height = "300px")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| MutualistBoth | Zr (mu) | 0.704 | 0.345 | 1.063 | 0.062 | 0.181 |
| MutualistHorizontal | Zr (mu) | 0.456 | 0.297 | 0.615 | NA | NA |
| MutualistVertical | Zr (mu) | 0.746 | 0.627 | 0.865 | NA | NA |
| ParasiteBoth | Zr (mu) | 0.564 | 0.442 | 0.686 | NA | NA |
| ParasiteHorizontal | Zr (mu) | 0.445 | 0.370 | 0.519 | NA | NA |
| MutualistBoth | r (mu) | 0.607 | 0.332 | 0.787 | NA | NA |
| MutualistHorizontal | r (mu) | 0.427 | 0.289 | 0.548 | NA | NA |
| MutualistVertical | r (mu) | 0.633 | 0.556 | 0.699 | NA | NA |
| ParasiteBoth | r (mu) | 0.511 | 0.415 | 0.595 | NA | NA |
| ParasiteHorizontal | r (mu) | 0.418 | 0.354 | 0.477 | NA | NA |
| MutualistBoth-MutualistHorizontal | Zr (beta) | -0.248 | -0.640 | 0.145 | NA | NA |
| MutualistBoth-MutualistVertical | Zr (beta) | 0.042 | -0.336 | 0.420 | NA | NA |
| MutualistBoth-ParasiteBoth | Zr (beta) | -0.140 | -0.519 | 0.239 | NA | NA |
| MutualistBoth-ParasiteHorizontal | Zr (beta) | -0.259 | -0.625 | 0.107 | NA | NA |
| MutualistHorizontal-MutualistVertical | Zr (beta) | 0.290 | 0.091 | 0.488 | NA | NA |
| MutualistHorizontal-ParasiteBoth | Zr (beta) | 0.108 | -0.093 | 0.308 | NA | NA |
| MutualistHorizontal-ParasiteHorizontal | Zr (beta) | -0.011 | -0.183 | 0.161 | NA | NA |
| MutualistVertical-ParasiteBoth | Zr (beta) | -0.182 | -0.352 | -0.012 | NA | NA |
| MutualistVertical-ParasiteHorizontal | Zr (beta) | -0.301 | -0.441 | -0.161 | NA | NA |
| ParasiteBoth-ParasiteHorizontal | Zr (beta) | -0.119 | -0.262 | 0.024 | NA | NA |
# colour list
colour_ls <- c("#000000", "#E69F00", "#56B4E9", "#009E73", "#F0E422", "#0072B2", "#D55E00", "#CC79A7", "#00008B", "#8B0A50", "#54FF9F", "#999999")
# adding sample size (k) for each category
k_symbiosis_transmission <- dat %>% group_by(symbiosis_transmission) %>% count()
# getting estimates and predicitons
pred_symbiosis_transmission <- get_pred(mr_symbiosis_transmission1, mod = "symbiosis_transmission")
res_symbiosis_transmission1 <- left_join(res_symbiosis_transmission1, k_symbiosis_transmission, by = c("name" = "symbiosis_transmission")) %>% left_join(pred_symbiosis_transmission)
#res_symbiosis1
# drawing a funnel plot - fig 2b
fig_symbiosis_transmission <- ggplot(data = res_symbiosis_transmission1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(symbiosis_transmission)),
aes(x= tanh(Zr), y = symbiosis_transmission, size = ((1/VZr) + 3), colour = symbiosis_transmission), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("MutualistBoth"= colour_ls[1], "MutualistHorizontal"= colour_ls[2], "MutualistVertical" = colour_ls[3],"ParasiteBoth"= colour_ls[4], "ParasiteHorizontal" = colour_ls[5])) +
scale_fill_manual(values = c("MutualistBoth"= colour_ls[1], "MutualistHorizontal"= colour_ls[2], "MutualistVertical" = colour_ls[3],"ParasiteBoth"= colour_ls[4], "ParasiteHorizontal" = colour_ls[5])) +
scale_y_discrete(labels = c("MutualistBoth" = "Mutualist-\nBoth", "MutualistHorizontal" = "Mutualist-\nHorizontal","MutualistVertical" = "Mutualist-\nVertical", "ParasiteBoth" = "Parasite-\nBoth", "ParasiteHorizontal" = "Parasite-\nHorizontal")) +
annotate('text', x = 0.93, y = 1:5 + 0.15, label= paste("italic(k)==", res_symbiosis_transmission1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")) ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position= c(0, 1), legend.justification = c(0,1)) +
theme(legend.direction="horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_mutualism), xmin = -1.1, xmax = -0.9, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_both), xmin = -0.9, xmax = -0.6, ymin = 0.4, ymax = 1.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -1.1, xmax = -0.9, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_horizontal),xmin = -0.9, xmax = -0.6, ymin = 1.4, ymax = 2.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -1.1, xmax = -0.9, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_vertical), xmin = -0.9, xmax = -0.6, ymin = 2.4, ymax = 3.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -1.1, xmax = -0.9, ymin = 3.6, ymax = 4.2) +
annotation_custom(rasterGrob(image_both), xmin = -0.9, xmax = -0.6, ymin = 3.4, ymax = 4.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -1.1, xmax = -0.9, ymin = 4.6, ymax = 5.2) +
annotation_custom(rasterGrob(image_horizontal), xmin = -0.9, xmax = -0.6, ymin = 4.4, ymax = 5.2)
fig_symbiosis_transmission
## fig 3
g <- ggplot(data = res_symbiosis_transmission1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(symbiosis_transmission)),
aes(x= tanh(Zr), y = symbiosis_transmission, size = ((1/VZr) + 3), colour = symbiosis_transmission), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("MutualistBoth"= colour_ls[1], "MutualistHorizontal"= colour_ls[2], "MutualistVertical" = colour_ls[3],"ParasiteBoth"= colour_ls[4], "ParasiteHorizontal" = colour_ls[5])) +
scale_fill_manual(values = c("MutualistBoth"= colour_ls[1], "MutualistHorizontal"= colour_ls[2], "MutualistVertical" = colour_ls[3],"ParasiteBoth"= colour_ls[4], "ParasiteHorizontal" = colour_ls[5])) +
scale_y_discrete(labels = c("MutualistBoth" = "Mutualist-\nBoth", "MutualistHorizontal" = "Mutualist-\nHorizontal","MutualistVertical" = "Mutualist-\nVertical", "ParasiteBoth" = "Parasite-\nBoth", "ParasiteHorizontal" = "Parasite-\nHorizontal")) +
annotate('text', x = 0.93, y = 1:5 + 0.15, label= paste("italic(k)==", res_symbiosis_transmission1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")) ,tag = "g" ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position="none") +
theme(axis.text.y = element_text(size = 10, colour ="black",hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_mutualism), xmin = -1.1, xmax = -0.9, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_both), xmin = -0.9, xmax = -0.6, ymin = 0.4, ymax = 1.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -1.1, xmax = -0.9, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_horizontal),xmin = -0.9, xmax = -0.6, ymin = 1.4, ymax = 2.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -1.1, xmax = -0.9, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_vertical), xmin = -0.9, xmax = -0.6, ymin = 2.4, ymax = 3.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -1.1, xmax = -0.9, ymin = 3.6, ymax = 4.2) +
annotation_custom(rasterGrob(image_both), xmin = -0.9, xmax = -0.6, ymin = 3.4, ymax = 4.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -1.1, xmax = -0.9, ymin = 4.6, ymax = 5.2) +
annotation_custom(rasterGrob(image_horizontal), xmin = -0.9, xmax = -0.6, ymin = 4.4, ymax = 5.2)
Figure 4c: An orchard plot showing the group-wise means (the categorical variable symbiosis_transmission) with their 95% confidence intervals (thick lines) and 95% prediction intervals (thin lines), with observed effect sizes based on various sample sizes.
# building fig 3 using patchwork
fig3 <- (a/b/c/d + plot_layout(heights = c(1.6, 2, 3.7, 3.7)))
# fig3 <- (a / b / c / d + plot_layout(heights = c(1.6, 2, 3.7, 3.7))) | (e / f /
# g + plot_layout(heights = c(2.8, 2.8, 4.4))) #+ plot_annotation(tag_levels =
# 'a', tag_suffix = ')')
fig3
# ggsave(here('figs/fig3.png'), width = 8, height = 12)
# ggsave(here('figs/fig3.pdf'), width = 8, height = 12)
Figure 3: putting all 4 panels together: Figure 3a - Figure 3d (see the main text)
# building fig 3 using patchwork
fig4 <- (e/f/g + plot_layout(heights = c(2.8, 2.8, 4.4))) + plot_annotation(tag_levels = "a")
fig4
# ggsave(here('figs/fig4.png'), width = 8, height = 12)
# ggsave(here('figs/fig4.pdf'), width = 8, height = 12)
Figure 4: putting all 3 panels together: Figure 4a - Figure 4c (see the main text)
These are extra analyses not discussed in the main text.
# reordering
dat$host_tax_symbiosis <- factor(dat$host_tax_symbiosis, levels = c("MicrobeMutualist",
"MicrobeParasite", "PlantMutualist", "PlantParasite", "InvertMutualist", "InvertParasite",
"VertMutualist", "VertParasite"))
# meta-regression: mutiple intercepts
mr_host_tax_symbiosis1 <- rma.mv(yi = Zr, V = VZr, mods = ~host_tax_symbiosis - 1,
test = "t", random = ~1 | authors, data = dat)
# # meta-regression: contrasts x 10 getting the level names out
level_names <- levels(dat$host_tax_symbiosis)
# helper function to run metafor meta-regression
run_rma <- function(name) {
rma.mv(yi = Zr, V = VZr, mods = ~relevel(host_tax_symbiosis, ref = name), test = "t",
random = ~1 | authors, data = dat)
}
# results of meta-regression including all contrast results; taking the last
# level out ([-length(level_names)])
mr_host_tax_symbiosis <- map(level_names[-length(level_names)], run_rma)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with host_tax_symbiosis. Note that mu means the group mean while beta represents the contrast between two groups in the Unit column.
# getting marginal R2
r2_host_tax_symbiosis1 <- R2(mr_host_tax_symbiosis1)
# getting estimates
res_host_tax_symbiosis1 <- get_est(mr_host_tax_symbiosis1, mod = "host_tax_symbiosis")
res_host_tax_symbiosis <- map(mr_host_tax_symbiosis, ~get_est(.x, mod = "host_tax_symbiosis"))
# a list of the numbers to take out unnecessary contrasts
contra_list <- Map(seq, from = 1, to = 1:7)
# you need to flatten twice: first to make it a list and make it a vector
estimates <- map2(res_host_tax_symbiosis, contra_list, ~.x[-(.y), "estimate"]) %>%
flatten() %>% flatten_dbl()
lowerCLs <- map2(res_host_tax_symbiosis, contra_list, ~.x[-(.y), "lowerCL"]) %>%
flatten() %>% flatten_dbl()
upperCLs <- map2(res_host_tax_symbiosis, contra_list, ~.x[-(.y), "upperCL"]) %>%
flatten() %>% flatten_dbl()
# creating a table
tibble(`Fixed effect` = c(as.character(res_host_tax_symbiosis1$name), cont_gen(res_host_tax_symbiosis1$name)),
Estimate = c(res_host_tax_symbiosis1$estimate, estimates), `Lower CI [0.025]` = c(res_host_tax_symbiosis1$lowerCL,
lowerCLs), `Upper CI [0.975]` = c(res_host_tax_symbiosis1$upperCL, upperCLs),
`V[authors]` = c(mr_host_tax_symbiosis1$sigma2, rep(NA, (8 + choose(8, 2)) -
1)), R2 = c(r2_host_tax_symbiosis1[1], rep(NA, (8 + choose(8, 2)) - 1))) %>%
kable("html", digits = 3) %>% kable_styling("striped", position = "left") %>%
scroll_box(width = "100%", height = "300px")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| MicrobeMutualist | 1.316 | 0.982 | 1.650 | 0.069 | 0.311 |
| MicrobeParasite | 0.419 | 0.082 | 0.755 | NA | NA |
| PlantMutualist | 0.373 | 0.197 | 0.548 | NA | NA |
| PlantParasite | 0.379 | 0.237 | 0.521 | NA | NA |
| InvertMutualist | 0.638 | 0.520 | 0.757 | NA | NA |
| InvertParasite | 0.552 | 0.336 | 0.768 | NA | NA |
| VertMutualist | 1.166 | 0.294 | 2.038 | NA | NA |
| VertParasite | 0.505 | 0.427 | 0.583 | NA | NA |
| MicrobeMutualist-MicrobeParasite | -0.898 | -1.372 | -0.423 | NA | NA |
| MicrobeMutualist-PlantMutualist | -0.943 | -1.320 | -0.566 | NA | NA |
| MicrobeMutualist-PlantParasite | -0.937 | -1.300 | -0.575 | NA | NA |
| MicrobeMutualist-InvertMutualist | -0.678 | -1.032 | -0.324 | NA | NA |
| MicrobeMutualist-InvertParasite | -0.764 | -1.162 | -0.367 | NA | NA |
| MicrobeMutualist-VertMutualist | -0.150 | -1.084 | 0.783 | NA | NA |
| MicrobeMutualist-VertParasite | -0.811 | -1.154 | -0.469 | NA | NA |
| MicrobeParasite-PlantMutualist | -0.046 | -0.425 | 0.334 | NA | NA |
| MicrobeParasite-PlantParasite | -0.040 | -0.405 | 0.326 | NA | NA |
| MicrobeParasite-InvertMutualist | 0.220 | -0.137 | 0.577 | NA | NA |
| MicrobeParasite-InvertParasite | 0.133 | -0.267 | 0.533 | NA | NA |
| MicrobeParasite-VertMutualist | 0.748 | -0.187 | 1.682 | NA | NA |
| MicrobeParasite-VertParasite | 0.086 | -0.259 | 0.432 | NA | NA |
| PlantMutualist-PlantParasite | 0.006 | -0.206 | 0.218 | NA | NA |
| PlantMutualist-InvertMutualist | 0.266 | 0.054 | 0.477 | NA | NA |
| PlantMutualist-InvertParasite | 0.179 | -0.099 | 0.457 | NA | NA |
| PlantMutualist-VertMutualist | 0.793 | -0.096 | 1.683 | NA | NA |
| PlantMutualist-VertParasite | 0.132 | -0.060 | 0.324 | NA | NA |
| PlantParasite-InvertMutualist | 0.260 | 0.075 | 0.444 | NA | NA |
| PlantParasite-InvertParasite | 0.173 | -0.085 | 0.431 | NA | NA |
| PlantParasite-VertMutualist | 0.787 | -0.096 | 1.671 | NA | NA |
| PlantParasite-VertParasite | 0.126 | -0.036 | 0.288 | NA | NA |
| InvertMutualist-InvertParasite | -0.087 | -0.333 | 0.159 | NA | NA |
| InvertMutualist-VertMutualist | 0.528 | -0.352 | 1.408 | NA | NA |
| InvertMutualist-VertParasite | -0.133 | -0.275 | 0.008 | NA | NA |
| InvertParasite-VertMutualist | 0.614 | -0.284 | 1.512 | NA | NA |
| InvertParasite-VertParasite | -0.047 | -0.270 | 0.176 | NA | NA |
| VertMutualist-VertParasite | -0.661 | -1.536 | 0.214 | NA | NA |
# adding sample size (k) for each category
k_host_tax_symbiosis <- dat %>% group_by(host_tax_symbiosis) %>% count()
# getting estimates and predicitons
pred_host_tax_symbiosis <- get_pred(mr_host_tax_symbiosis1, mod = "host_tax_symbiosis")
res_host_tax_symbiosis1 <- left_join(res_host_tax_symbiosis1, k_host_tax_symbiosis, by = c("name" = "host_tax_symbiosis")) %>% left_join(pred_host_tax_symbiosis)
#res_symbiosis1
# drawing a funnel plot - fig 2b
fig_host_tax_symbiosis <- ggplot(data = res_host_tax_symbiosis1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(host_tax_symbiosis)),
aes(x= tanh(Zr), y = host_tax_symbiosis, size = ((1/VZr) + 3), colour = host_tax_symbiosis), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("MicrobeMutualist"= colour_ls[1], "MicrobeParasite"= colour_ls[2], "PlantMutualist"= colour_ls[3], "PlantParasite"= colour_ls[4], "InvertMutualist" = colour_ls[5], "InvertParasite"= colour_ls[6], "VertMutualist"= colour_ls[7], "VertParasite"= colour_ls[8] )) +
scale_fill_manual(values = c("MicrobeMutualist"= colour_ls[1], "MicrobeParasite"= colour_ls[2], "PlantMutualist"= colour_ls[3], "PlantParasite"= colour_ls[4], "InvertMutualist" = colour_ls[5], "InvertParasite"= colour_ls[6], "VertMutualist"= colour_ls[7], "VertParasite"= colour_ls[8] )) +
scale_y_discrete(labels = c("MicrobeMutualist"= "Microbe-\nMutualist", "MicrobeParasite"= "Microbe-\nParasite", "PlantMutualist" = "Plant-\nMutualist", "PlantParasite"="Plant-\nParasite", "InvertMutualist" = "Invertebrate-\nMutualist", "InvertParasite"= "Invertebrate-\nParasite", "VertMutualist"= "Vertebrate-\nMutualist", "VertParasite"= "Vertebrate-\nParasite" )) +
annotate('text', x = 0.93, y = 1:8 + 0.15, label= paste("italic(k)==", res_host_tax_symbiosis1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")) ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position= c(0, 1), legend.justification = c(0,1)) +
theme(legend.direction="horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_microbe_host), xmin = -1.1, xmax = -0.9, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -0.9, xmax = -0.7, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_microbe_host), xmin = -1.1, xmax = -0.9, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_parasitism),xmin = -0.9, xmax = -0.7, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_plant_host), xmin = -1.1, xmax = -0.9, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -0.9, xmax = -0.7, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_plant_host), xmin = -1.1, xmax = -0.9, ymin = 3.6, ymax = 4.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -0.9, xmax = -0.7, ymin = 3.6, ymax = 4.2) +
annotation_custom(rasterGrob(image_invertebrate_host), xmin = -1.1, xmax = -0.9, ymin = 4.6, ymax = 5.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -0.9, xmax = -0.7, ymin = 4.6, ymax = 5.2) +
annotation_custom(rasterGrob(image_invertebrate_host), xmin = -1.1, xmax = -0.9, ymin = 5.6, ymax = 6.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -0.9, xmax = -0.7, ymin = 5.6, ymax = 6.2) +
annotation_custom(rasterGrob(image_vertebrate_host), xmin = -1.1, xmax = -0.9, ymin = 6.6, ymax = 7.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -0.9, xmax = -0.7, ymin = 6.6, ymax = 7.2) +
annotation_custom(rasterGrob(image_vertebrate_host), xmin = -1.1, xmax = -0.9, ymin = 7.6, ymax = 8.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -0.9, xmax = -0.7, ymin = 7.6, ymax = 8.2)
fig_host_tax_symbiosis
Figure S3.3: An orchard plot showing group-wise means (the categorical variable host_tax_symbiosis) with their 95% confidence intervals (thick lines) and 95% prediction intervals (thin lines), with observed effect sizes based on various sample sizes.
Splitting host taxonomy by mode of symbiosis revealed that the observed higher phylogenetic congruence of host-symbiont cophylogenies involving a microbial host is driven primarily by greater congruence between microbial hosts and mutualist symbionts. Congruence is also relatively high for invertebrate hosts that harbour a mutualistic symbiont, while congruence appears to be lowest for plant hosts that harbour a parasitic symbiont.
# reordering
dat$symbiont_tax_symbiosis <- factor(dat$symbiont_tax_symbiosis, levels = c("MicrobeMutualist",
"MicrobeParasite", "PlantMutualist", "PlantParasite", "InvertMutualist", "InvertParasite",
"VertParasite"))
# meta-regression: multiple intercepts
mr_symbiont_tax_symbiosis1 <- rma.mv(yi = Zr, V = VZr, mods = ~symbiont_tax_symbiosis -
1, test = "t", random = ~1 | authors, data = dat)
# # meta-regression: contrasts x 10 getting the level names out
level_names <- levels(dat$symbiont_tax_symbiosis)
# helper function to run metafor meta-regression
run_rma <- function(name) {
rma.mv(yi = Zr, V = VZr, mods = ~relevel(symbiont_tax_symbiosis, ref = name),
test = "t", random = ~1 | authors, data = dat)
}
# results of meta-regression including all contrast results; taking the last
# level out ([-length(level_names)])
mr_symbiont_tax_symbiosis <- map(level_names[-length(level_names)], run_rma)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with symbiont_tax_symbiosis. Note that mu means the group mean while beta represents the contrast between two groups in the Unit column.
# getting marginal R2
r2_symbiont_tax_symbiosis1 <- R2(mr_symbiont_tax_symbiosis1)
# getting estimates
res_symbiont_tax_symbiosis1 <- get_est(mr_symbiont_tax_symbiosis1, mod = "symbiont_tax_symbiosis")
res_symbiont_tax_symbiosis <- map(mr_symbiont_tax_symbiosis, ~get_est(.x, mod = "symbiont_tax_symbiosis"))
# a list of the numbers to take out unnecessary contrasts
contra_list <- Map(seq, from = 1, to = 1:6)
# you need to flatten twice: first to make it a list and make it a vector
estimates <- map2(res_symbiont_tax_symbiosis, contra_list, ~.x[-(.y), "estimate"]) %>%
flatten() %>% flatten_dbl()
lowerCLs <- map2(res_symbiont_tax_symbiosis, contra_list, ~.x[-(.y), "lowerCL"]) %>%
flatten() %>% flatten_dbl()
upperCLs <- map2(res_symbiont_tax_symbiosis, contra_list, ~.x[-(.y), "upperCL"]) %>%
flatten() %>% flatten_dbl()
# creating a table
tibble(`Fixed effect` = c(as.character(res_symbiont_tax_symbiosis1$name), cont_gen(res_symbiont_tax_symbiosis1$name)),
Estimate = c(res_symbiont_tax_symbiosis1$estimate, estimates), `Lower CI [0.025]` = c(res_symbiont_tax_symbiosis1$lowerCL,
lowerCLs), `Upper CI [0.975]` = c(res_symbiont_tax_symbiosis1$upperCL, upperCLs),
`V[authors]` = c(mr_symbiont_tax_symbiosis1$sigma2, rep(NA, (7 + choose(7, 2)) -
1)), R2 = c(r2_symbiont_tax_symbiosis1[1], rep(NA, (7 + choose(7, 2)) - 1))) %>%
kable("html", digits = 3) %>% kable_styling("striped", position = "left") %>%
scroll_box(width = "100%", height = "300px")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| MicrobeMutualist | 0.670 | 0.557 | 0.782 | 0.075 | 0.389 |
| MicrobeParasite | 0.435 | 0.330 | 0.539 | NA | NA |
| PlantMutualist | 1.009 | 0.536 | 1.482 | NA | NA |
| PlantParasite | 3.453 | 1.411 | 5.496 | NA | NA |
| InvertMutualist | 0.445 | 0.248 | 0.641 | NA | NA |
| InvertParasite | 0.510 | 0.422 | 0.597 | NA | NA |
| VertParasite | 0.440 | -0.173 | 1.053 | NA | NA |
| MicrobeMutualist-MicrobeParasite | -0.235 | -0.388 | -0.082 | NA | NA |
| MicrobeMutualist-PlantMutualist | 0.339 | -0.147 | 0.826 | NA | NA |
| MicrobeMutualist-PlantParasite | 2.784 | 0.738 | 4.830 | NA | NA |
| MicrobeMutualist-InvertMutualist | -0.225 | -0.452 | 0.002 | NA | NA |
| MicrobeMutualist-InvertParasite | -0.160 | -0.303 | -0.018 | NA | NA |
| MicrobeMutualist-VertParasite | -0.230 | -0.853 | 0.394 | NA | NA |
| MicrobeParasite-PlantMutualist | 0.574 | 0.090 | 1.059 | NA | NA |
| MicrobeParasite-PlantParasite | 3.019 | 0.973 | 5.064 | NA | NA |
| MicrobeParasite-InvertMutualist | 0.010 | -0.213 | 0.233 | NA | NA |
| MicrobeParasite-InvertParasite | 0.075 | -0.061 | 0.211 | NA | NA |
| MicrobeParasite-VertParasite | 0.005 | -0.616 | 0.627 | NA | NA |
| PlantMutualist-PlantParasite | 2.444 | 0.347 | 4.541 | NA | NA |
| PlantMutualist-InvertMutualist | -0.564 | -1.077 | -0.052 | NA | NA |
| PlantMutualist-InvertParasite | -0.499 | -0.981 | -0.018 | NA | NA |
| PlantMutualist-VertParasite | -0.569 | -1.343 | 0.205 | NA | NA |
| PlantParasite-InvertMutualist | -3.009 | -5.061 | -0.956 | NA | NA |
| PlantParasite-InvertParasite | -2.944 | -4.989 | -0.899 | NA | NA |
| PlantParasite-VertParasite | -3.013 | -5.146 | -0.880 | NA | NA |
| InvertMutualist-InvertParasite | 0.065 | -0.143 | 0.273 | NA | NA |
| InvertMutualist-VertParasite | -0.005 | -0.648 | 0.639 | NA | NA |
| InvertParasite-VertParasite | -0.070 | -0.689 | 0.550 | NA | NA |
# adding sample size (k) for each category
k_symbiont_tax_symbiosis <- dat %>% group_by(symbiont_tax_symbiosis) %>% count()
# getting estimates and predicitons
pred_symbiont_tax_symbiosis <- get_pred(mr_symbiont_tax_symbiosis1, mod = "symbiont_tax_symbiosis")
res_symbiont_tax_symbiosis1 <- left_join(res_symbiont_tax_symbiosis1, k_symbiont_tax_symbiosis, by = c("name" = "symbiont_tax_symbiosis")) %>% left_join(pred_symbiont_tax_symbiosis)
#res_symbiosis1
# drawing a funnel plot - fig 2b
fig_symbiont_tax_symbiosis <- ggplot(data = res_symbiont_tax_symbiosis1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(symbiont_tax_symbiosis)),
aes(x= tanh(Zr), y = symbiont_tax_symbiosis, size = ((1/VZr) + 3), colour = symbiont_tax_symbiosis), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
# setting colours
scale_color_manual(values = c("MicrobeMutualist"= colour_ls[1], "MicrobeParasite"= colour_ls[2], "PlantMutualist"= colour_ls[3], "PlantParasite"= colour_ls[4], "InvertMutualist" = colour_ls[5], "InvertParasite"= colour_ls[6], "VertParasite"= colour_ls[8] )) +
scale_fill_manual(values = c("MicrobeMutualist"= colour_ls[1], "MicrobeParasite"= colour_ls[2], "PlantMutualist"= colour_ls[3], "PlantParasite"= colour_ls[4], "InvertMutualist" = colour_ls[5], "InvertParasite"= colour_ls[6], "VertParasite"= colour_ls[8] )) +
scale_y_discrete(labels = c("MicrobeMutualist"= "Microbe-\nMutualist", "MicrobeParasite"= "Microbe-\nParasite", "PlantMutualist" = "Plant-\nMutualist", "PlantParasite"="Plant-\nParasite", "InvertMutualist" = "Invertebrate-\nMutualist", "InvertParasite"= "Invertebrate-\nParasite", "VertParasite"= "Vertebrate-\nParasite" )) +
annotate('text', x = 0.93, y = 1:7 + 0.15, label= paste("italic(k)==", res_symbiont_tax_symbiosis1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")) ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position= c(0, 1), legend.justification = c(0, 1)) +
theme(legend.direction="horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_microbe_parasite), xmin = -1.1, xmax = -0.9, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -0.9, xmax = -0.7, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_microbe_parasite), xmin = -1.1, xmax = -0.9, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_parasitism),xmin = -0.9, xmax = -0.7, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_plant_parasite), xmin = -1.1, xmax = -0.9, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -0.9, xmax = -0.7, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_plant_parasite), xmin = -1.1, xmax = -0.9, ymin = 3.6, ymax = 4.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -0.9, xmax = -0.7, ymin = 3.6, ymax = 4.2) +
annotation_custom(rasterGrob(image_invertebrate_parasite), xmin = -1.1, xmax = -0.9, ymin = 4.6, ymax = 5.2) +
annotation_custom(rasterGrob(image_mutualism), xmin = -0.9, xmax = -0.7, ymin = 4.6, ymax = 5.2) +
annotation_custom(rasterGrob(image_invertebrate_parasite), xmin = -1.1, xmax = -0.9, ymin = 5.6, ymax = 6.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -0.9, xmax = -0.7, ymin = 5.6, ymax = 6.2) +
annotation_custom(rasterGrob(image_vertebrate_parasite), xmin = -1.1, xmax = -0.9, ymin = 6.6, ymax = 7.2) +
annotation_custom(rasterGrob(image_parasitism), xmin = -0.9, xmax = -0.7, ymin = 6.6, ymax = 7.2)
fig_symbiont_tax_symbiosis
Figure S3.4: An orchard plot showing group-wise means (the categorical variable symbiont_tax_symbiosis) with their 95% confidence intervals (thick lines) and 95% prediction intervals (thin lines), with observed effect sizes based on various sample sizes.
Splitting symbiont taxonomy by mode of symbiosis revealed much less variation, except for higher congruence exhibited by cophylogenies involving a plant symbiont (instances of which are relatively rare), and the finding that cophylogenies involving a microbial mutualist symbiont are slightly more congruent than the remaining categories.
# reordering
dat$host_symbiont_tax <- factor(dat$host_symbiont_tax, levels = c("MicrobeInvert",
"MicrobeMicrobe", "MicrobePlant", "PlantInvert", "PlantMicrobe", "InvertInvert",
"InvertMicrobe", "InvertPlant", "VertInvert", "VertMicrobe", "VertVert"))
# meta-regression: multiple intercepts
mr_host_symbiont_tax1 <- rma.mv(yi = Zr, V = VZr, mods = ~host_symbiont_tax - 1,
test = "t", random = ~1 | authors, data = dat)
# # meta-regression: contrasts x 10 getting the level names out
level_names <- levels(dat$host_symbiont_tax)
# helper function to run metafor meta-regression
run_rma <- function(name) {
rma.mv(yi = Zr, V = VZr, mods = ~relevel(host_symbiont_tax, ref = name), test = "t",
random = ~1 | authors, data = dat)
}
# results of meta-regression including all contrast results; taking the last
# level out ([-length(level_names)])
mr_host_symbiont_tax <- map(level_names[-length(level_names)], run_rma)
Regression coefficients (estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with host_symbiont_tax. Note that mu means the group mean while beta represents the contrast between two groups in the Unit column.
# getting marginal R2
r2_host_symbiont_tax1 <- R2(mr_host_symbiont_tax1)
# getting estimates
res_host_symbiont_tax1 <- get_est(mr_host_symbiont_tax1, mod = "host_symbiont_tax")
res_host_symbiont_tax <- map(mr_host_symbiont_tax, ~get_est(.x, mod = "host_symbiont_tax"))
# a list of the numbers to take out unnecessary contrasts
contra_list <- Map(seq, from = 1, to = 1:10)
# you need to flatten twice: first to make it a list and make it a vector
estimates <- map2(res_host_symbiont_tax, contra_list, ~.x[-(.y), "estimate"]) %>%
flatten() %>% flatten_dbl()
lowerCLs <- map2(res_host_symbiont_tax, contra_list, ~.x[-(.y), "lowerCL"]) %>% flatten() %>%
flatten_dbl()
upperCLs <- map2(res_host_symbiont_tax, contra_list, ~.x[-(.y), "upperCL"]) %>% flatten() %>%
flatten_dbl()
# creating a table
tibble(`Fixed effect` = c(as.character(res_host_symbiont_tax1$name), cont_gen(res_host_symbiont_tax1$name)),
Estimate = c(res_host_symbiont_tax1$estimate, estimates), `Lower CI [0.025]` = c(res_host_symbiont_tax1$lowerCL,
lowerCLs), `Upper CI [0.975]` = c(res_host_symbiont_tax1$upperCL, upperCLs),
`V[authors]` = c(mr_host_tax_symbiosis1$sigma2, rep(NA, (11 + choose(11, 2)) -
1)), R2 = c(r2_host_symbiont_tax1[1], rep(NA, (11 + choose(11, 2)) - 1))) %>%
kable("html", digits = 3) %>% kable_styling("striped", position = "left") %>%
scroll_box(width = "100%", height = "300px")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| MicrobeInvert | 1.069 | -0.196 | 2.335 | 0.069 | 0.23 |
| MicrobeMicrobe | 0.820 | 0.532 | 1.107 | NA | NA |
| MicrobePlant | 1.055 | 0.552 | 1.558 | NA | NA |
| PlantInvert | 0.329 | 0.176 | 0.482 | NA | NA |
| PlantMicrobe | 0.460 | 0.258 | 0.662 | NA | NA |
| InvertInvert | 0.661 | 0.449 | 0.874 | NA | NA |
| InvertMicrobe | 0.612 | 0.487 | 0.738 | NA | NA |
| InvertPlant | 1.665 | 0.399 | 2.931 | NA | NA |
| VertInvert | 0.562 | 0.454 | 0.669 | NA | NA |
| VertMicrobe | 0.459 | 0.335 | 0.584 | NA | NA |
| VertVert | 0.440 | -0.187 | 1.067 | NA | NA |
| MicrobeInvert-MicrobeMicrobe | -0.250 | -1.547 | 1.048 | NA | NA |
| MicrobeInvert-MicrobePlant | -0.014 | -1.376 | 1.348 | NA | NA |
| MicrobeInvert-PlantInvert | -0.740 | -2.015 | 0.535 | NA | NA |
| MicrobeInvert-PlantMicrobe | -0.609 | -1.891 | 0.672 | NA | NA |
| MicrobeInvert-InvertInvert | -0.408 | -1.691 | 0.875 | NA | NA |
| MicrobeInvert-InvertMicrobe | -0.457 | -1.728 | 0.815 | NA | NA |
| MicrobeInvert-InvertPlant | 0.596 | -1.194 | 2.386 | NA | NA |
| MicrobeInvert-VertInvert | -0.507 | -1.777 | 0.763 | NA | NA |
| MicrobeInvert-VertMicrobe | -0.610 | -1.881 | 0.662 | NA | NA |
| MicrobeInvert-VertVert | -0.629 | -2.041 | 0.783 | NA | NA |
| MicrobeMicrobe-MicrobePlant | 0.236 | -0.343 | 0.815 | NA | NA |
| MicrobeMicrobe-PlantInvert | -0.491 | -0.816 | -0.165 | NA | NA |
| MicrobeMicrobe-PlantMicrobe | -0.360 | -0.711 | -0.008 | NA | NA |
| MicrobeMicrobe-InvertInvert | -0.158 | -0.516 | 0.199 | NA | NA |
| MicrobeMicrobe-InvertMicrobe | -0.207 | -0.520 | 0.106 | NA | NA |
| MicrobeMicrobe-InvertPlant | 0.845 | -0.452 | 2.143 | NA | NA |
| MicrobeMicrobe-VertInvert | -0.258 | -0.564 | 0.049 | NA | NA |
| MicrobeMicrobe-VertMicrobe | -0.360 | -0.673 | -0.047 | NA | NA |
| MicrobeMicrobe-VertVert | -0.379 | -1.069 | 0.310 | NA | NA |
| MicrobePlant-PlantInvert | -0.726 | -1.252 | -0.201 | NA | NA |
| MicrobePlant-PlantMicrobe | -0.595 | -1.137 | -0.053 | NA | NA |
| MicrobePlant-InvertInvert | -0.394 | -0.940 | 0.152 | NA | NA |
| MicrobePlant-InvertMicrobe | -0.443 | -0.961 | 0.075 | NA | NA |
| MicrobePlant-InvertPlant | 0.610 | -0.752 | 1.972 | NA | NA |
| MicrobePlant-VertInvert | -0.493 | -1.007 | 0.021 | NA | NA |
| MicrobePlant-VertMicrobe | -0.596 | -1.114 | -0.078 | NA | NA |
| MicrobePlant-VertVert | -0.615 | -1.418 | 0.188 | NA | NA |
| PlantInvert-PlantMicrobe | 0.131 | -0.122 | 0.384 | NA | NA |
| PlantInvert-InvertInvert | 0.332 | 0.070 | 0.594 | NA | NA |
| PlantInvert-InvertMicrobe | 0.283 | 0.086 | 0.481 | NA | NA |
| PlantInvert-InvertPlant | 1.336 | 0.061 | 2.611 | NA | NA |
| PlantInvert-VertInvert | 0.233 | 0.046 | 0.420 | NA | NA |
| PlantInvert-VertMicrobe | 0.130 | -0.067 | 0.327 | NA | NA |
| PlantInvert-VertVert | 0.111 | -0.534 | 0.756 | NA | NA |
| PlantMicrobe-InvertInvert | 0.201 | -0.092 | 0.495 | NA | NA |
| PlantMicrobe-InvertMicrobe | 0.152 | -0.085 | 0.390 | NA | NA |
| PlantMicrobe-InvertPlant | 1.205 | -0.077 | 2.487 | NA | NA |
| PlantMicrobe-VertInvert | 0.102 | -0.127 | 0.331 | NA | NA |
| PlantMicrobe-VertMicrobe | -0.001 | -0.238 | 0.237 | NA | NA |
| PlantMicrobe-VertVert | -0.020 | -0.678 | 0.638 | NA | NA |
| InvertInvert-InvertMicrobe | -0.049 | -0.296 | 0.198 | NA | NA |
| InvertInvert-InvertPlant | 1.004 | -0.280 | 2.287 | NA | NA |
| InvertInvert-VertInvert | -0.099 | -0.326 | 0.127 | NA | NA |
| InvertInvert-VertMicrobe | -0.202 | -0.448 | 0.044 | NA | NA |
| InvertInvert-VertVert | -0.221 | -0.883 | 0.441 | NA | NA |
| InvertMicrobe-InvertPlant | 1.053 | -0.219 | 2.324 | NA | NA |
| InvertMicrobe-VertInvert | -0.051 | -0.216 | 0.114 | NA | NA |
| InvertMicrobe-VertMicrobe | -0.153 | -0.329 | 0.023 | NA | NA |
| InvertMicrobe-VertVert | -0.172 | -0.811 | 0.467 | NA | NA |
| InvertPlant-VertInvert | -1.103 | -2.373 | 0.167 | NA | NA |
| InvertPlant-VertMicrobe | -1.206 | -2.477 | 0.066 | NA | NA |
| InvertPlant-VertVert | -1.225 | -2.637 | 0.187 | NA | NA |
| VertInvert-VertMicrobe | -0.103 | -0.267 | 0.062 | NA | NA |
| VertInvert-VertVert | -0.122 | -0.758 | 0.514 | NA | NA |
| VertMicrobe-VertVert | -0.019 | -0.658 | 0.620 | NA | NA |
# colour list
#colour_ls <- c("#000000", "#E69F00", "#56B4E9", "#009E73", "#F0E422", "#0072B2", "#D55E00", "#CC79A7", "#00008B", "#8B0A50", "#54FF9F", "#999999")
# adding sample size (k) for each category
k_host_symbiont_tax <- dat %>% group_by(host_symbiont_tax) %>% count()
# getting estimates and predicitons
pred_host_symbiont_tax <- get_pred(mr_host_symbiont_tax1, mod = "host_symbiont_tax")
res_host_symbiont_tax1 <- left_join(res_host_symbiont_tax1, k_host_symbiont_tax, by = c("name" = "host_symbiont_tax")) %>% left_join(pred_host_symbiont_tax)
#res_symbiosis1
# drawing a funnel plot - fig 2b
fig_host_symbiont_tax <- ggplot(data = res_host_symbiont_tax1, aes(x = tanh(estimate), y = name)) +
scale_x_continuous(limits=c(-1, 1), breaks = seq(-1, 1, by = 0.2) ) +
geom_quasirandom(data = dat %>% filter(!is.na(host_symbiont_tax)),
aes(x= tanh(Zr), y = host_symbiont_tax, size = ((1/VZr) + 3), colour = host_symbiont_tax), groupOnX = FALSE, alpha=0.4) +
# 95 %precition interval (PI)
geom_errorbarh(aes(xmin = tanh(lowerPR), xmax = tanh(upperPR)), height = 0, show.legend = F, size = 0.5, alpha = 0.6) +
# 95 %CI
geom_errorbarh(aes(xmin = tanh(lowerCL), xmax = tanh(upperCL)), height = 0, show.legend = F, size = 1.2) +
geom_vline(xintercept = 0, linetype = 2, colour = "black", alpha = 0.3) +
# creating dots and different size (bee-swarm and bubbles)
geom_point(aes(fill = name), size = 3, shape = 21) + #
# setting colours
scale_color_manual(values = c("MicrobeInvert" = colour_ls[1], "MicrobeMicrobe"= colour_ls[2], "MicrobePlant" = colour_ls[3], "PlantInvert" = colour_ls[4],"PlantMicrobe" = colour_ls[5], "InvertInvert" = colour_ls[6], "InvertMicrobe" = colour_ls[7], "InvertPlant" = colour_ls[8],"VertInvert" = colour_ls[9], "VertMicrobe"= colour_ls[10],"VertVert" = colour_ls[11])) +
scale_fill_manual(values = c("MicrobeInvert" = colour_ls[1], "MicrobeMicrobe"= colour_ls[2], "MicrobePlant" = colour_ls[3], "PlantInvert" = colour_ls[4],"PlantMicrobe" = colour_ls[5], "InvertInvert" = colour_ls[6], "InvertMicrobe" = colour_ls[7], "InvertPlant" = colour_ls[8],"VertInvert" = colour_ls[9], "VertMicrobe"= colour_ls[10],"VertVert" = colour_ls[11])) +
scale_y_discrete(labels = c("MicrobeInvert" = "Microbe-\nInvertebrate", "MicrobeMicrobe"= "Microbe-\nMicrobe", "MicrobePlant" = "Microbe-\nPlant", "PlantInvert" = "Plant-\nInvertebrate","PlantMicrobe" = "Plant-\nMicrobe", "InvertInvert" = "Invertebrate\nInvertebrate", "InvertMicrobe" = "Invertebrate-\nMicrobe", "InvertPlant" = "Invertebrate-\nPlant","VertInvert" = "Vertebrate-\nInvertebrate", "VertMicrobe"= "Vertebrate-\nMicrobe", "VertVert" = "Vertebrate-\nVertebrate")) +
annotate('text', x = 0.93, y = 1:11 + 0.15, label= paste("italic(k)==", res_host_symbiont_tax1$n), parse=TRUE, hjust = "left", size=3.5) +
labs(x = expression(paste(italic(r), " (correlation)")), y = "", size = expression(paste(italic(n), " (# of species pairs)")) ) +
guides(fill = "none", colour = "none") +
theme_bw() +
theme(legend.position= c(0, 1), legend.justification = c(0,1)) +
theme(legend.direction="horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90)) +
# putting pictures in
annotation_custom(rasterGrob(image_microbe_host), xmin = -1.1, xmax = -0.9, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_invertebrate_parasite), xmin = -0.9, xmax = -0.7, ymin = 0.6, ymax = 1.2) +
annotation_custom(rasterGrob(image_microbe_host), xmin = -1.1, xmax = -0.9, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_microbe_parasite),xmin = -0.9, xmax = -0.7, ymin = 1.6, ymax = 2.2) +
annotation_custom(rasterGrob(image_microbe_host), xmin = -1.1, xmax = -0.9, ymin = 2.6, ymax = 3.2) +
annotation_custom(rasterGrob(image_plant_parasite), xmin = -0.9, xmax = -0.7, ymin = 2.6, ymax = 3.2) +
#
annotation_custom(rasterGrob(image_plant_host), xmin = -1.1, xmax = -0.9, ymin = 3.6, ymax = 4.2) +
annotation_custom(rasterGrob(image_invertebrate_parasite), xmin = -0.9, xmax = -0.7, ymin = 3.6, ymax = 4.2) +
annotation_custom(rasterGrob(image_plant_host), xmin = -1.1, xmax = -0.9, ymin = 4.6, ymax = 5.2) +
annotation_custom(rasterGrob(image_microbe_parasite), xmin = -0.9, xmax = -0.7, ymin = 4.6, ymax = 5.2) +
#
annotation_custom(rasterGrob(image_invertebrate_host), xmin = -1.1, xmax = -0.9, ymin = 5.6, ymax = 6.2) +
annotation_custom(rasterGrob(image_invertebrate_parasite), xmin = -0.9, xmax = -0.7, ymin = 5.6, ymax = 6.2) +
annotation_custom(rasterGrob(image_invertebrate_host), xmin = -1.1, xmax = -0.9, ymin = 6.6, ymax = 7.2) +
annotation_custom(rasterGrob(image_microbe_parasite), xmin = -0.9, xmax = -0.7, ymin = 6.6, ymax = 7.2) +
annotation_custom(rasterGrob(image_invertebrate_host), xmin = -1.1, xmax = -0.9, ymin = 7.6, ymax = 8.2) +
annotation_custom(rasterGrob(image_plant_parasite), xmin = -0.9, xmax = -0.7, ymin = 7.6, ymax = 8.2) +
#
annotation_custom(rasterGrob(image_vertebrate_host), xmin = -1.1, xmax = -0.9, ymin = 8.6, ymax = 9.2) +
annotation_custom(rasterGrob(image_invertebrate_parasite), xmin = -0.9, xmax = -0.7, ymin = 8.6, ymax = 9.2) +
annotation_custom(rasterGrob(image_vertebrate_host), xmin = -1.1, xmax = -0.9, ymin = 9.6, ymax = 10.2) +
annotation_custom(rasterGrob(image_microbe_parasite), xmin = -0.9, xmax = -0.7, ymin = 9.6, ymax = 10.2) +
annotation_custom(rasterGrob(image_vertebrate_host), xmin = -1.1, xmax = -0.9, ymin = 10.6, ymax = 11.2) +
annotation_custom(rasterGrob(image_vertebrate_parasite), xmin = -0.9, xmax = -0.7, ymin = 10.6, ymax = 11.2)
fig_host_symbiont_tax
Figure S3.5: An orchard plot showing the group-wise means (the categorical variable host_symbiont_tax) with their 95% confidence intervals (thick lines) and 95% prediction intervals (thin lines), with observed effect sizes based on various sample sizes.
Here we build the best model via an AICc based model selection method implemented in the R package MuMin(Barton 2009). For the full model, we had 6 variables: symbiosis, host_tax_broad, symbiont_tax_broad, mode_of_transmission_broad, endo_or_ecto, & log(host_range_link_ratio). We did not use log(host_range_taxonomic_breadth) as it is co-linear with log(host_range_link_ratio) and also many of the interaction terms.
# creates a new function to run in MuMIn
updated.rma.mv <- updateable(rma.mv)
# updated.rma.mv
# testing the new function use method = 'ML' so that we can compare AIC
mr_full <- updated.rma.mv(yi = Zr, V = VZr, mods = ~symbiosis + host_tax_broad +
symbiont_tax_broad + mode_of_transmission_broad + endo_or_ecto + log(host_range_link_ratio),
test = "t", random = ~1 | authors, method = "ML", data = dat)
# ============================= additional methods for 'rma.mv' class (made by
# Kamil Barton) we need this to run model selection with rma.mv in MuMIn
# =============================
formula.rma.mv <- function(x, ...) return(eval(getCall(x)$mods))
makeArgs.rma.mv <- function(obj, termNames, comb, opt, ...) {
ret <- MuMIn:::makeArgs.default(obj, termNames, comb, opt)
names(ret)[1L] <- "mods"
ret
}
nobs.rma.mv <- function(object, ...) attr(logLik(object), "nall")
coefTable.rma.mv <- function(model, ...) MuMIn:::.makeCoefTable(model$b, model$se,
coefNames = rownames(model$b))
# =============================
# testing dredge dredge(full.model, evaluate=F) # show all candidate models n =
# 32 model exisit
candidates <- dredge(mr_full)
# displays delta AICc <2
candidates_aic2 <- subset(candidates, delta < 2)
# model averaging it seems like models are using z values rather than t values
# (which will be OK)
mr_averaged_aic2 <- summary(model.avg(candidates, delta < 2))
# relative importance of each predictor
importance <- importance(candidates)
# use REML if not for model comparision
model1 <- rma.mv(yi = Zr, V = VZr, mods = ~host_tax_broad + mode_of_transmission_broad +
symbiosis, test = "t", random = ~1 | authors, method = "REML", data = dat)
model2 <- rma.mv(yi = Zr, V = VZr, mods = ~host_tax_broad + mode_of_transmission_broad,
test = "t", random = ~1 | authors, method = "REML", data = dat)
model3 <- rma.mv(yi = Zr, V = VZr, mods = ~host_tax_broad + log(host_range_link_ratio) +
mode_of_transmission_broad + symbiosis, test = "t", random = ~1 | authors, method = "REML",
data = dat)
model4 <- rma.mv(yi = Zr, V = VZr, mods = ~host_tax_broad + log(host_range_link_ratio) +
mode_of_transmission_broad, test = "t", random = ~1 | authors, method = "REML",
data = dat)
The top 4 models (out of 32 possible models) within the \(\Delta\)AIC difference of 2, and which 6 variables: symbiosis, host_tax_broad, symbiont_tax_broad, mode_of_transmission_broad, endo_or_ecto, & log(host_range_link_ratio) were included (indicated by \(+\)); model weights (for the 2 models) and the sum of weights for each of the variables (from the 32 models) are included.
# creating a table
tibble(`Model (variable weight)` = c("Model1", "Model2", "Model3", "Model4", "(Sum of weights)"),
transmission = c(if_else(candidates_aic2$mode_of_transmission_broad == "+", "$+$",
"NA"), round(importance[1], 3)), host_tax = c(if_else(candidates_aic2$host_tax_broad ==
"+", "$+$", "NA"), round(importance[2], 3)), symbiosis = c(if_else(candidates_aic2$symbiosis ==
"+", "$+$", "NA"), round(importance[3], 3)), host_range = c(if_else(candidates_aic2$`log(host_range_link_ratio)` >=
0, "$+$", "NA"), round(importance[4], 3)), symbiont_tax = c(if_else(candidates_aic2$symbiont_tax_broad ==
"+", "$+$", "NA"), round(importance[5], 3)), endo_or_ecto = c(if_else(candidates_aic2$endo_or_ecto ==
"+", "$+$", "NA"), round(importance[6], 3)), delta_AICc = c(candidates_aic2$delta,
NA), Weight = c(candidates_aic2$weight, NA)) %>% kable("html", digits = 3) %>%
kable_styling("striped", position = "left")
| Model (variable weight) | transmission | host_tax | symbiosis | host_range | symbiont_tax | endo_or_ecto | delta_AICc | Weight |
|---|---|---|---|---|---|---|---|---|
| Model1 | \(+\) | \(+\) | \(+\) | NA | NA | NA | 0.000 | 0.369 |
| Model2 | \(+\) | \(+\) | NA | NA | NA | NA | 0.375 | 0.305 |
| Model3 | \(+\) | \(+\) | \(+\) | \(+\) | NA | NA | 1.520 | 0.172 |
| Model4 | \(+\) | \(+\) | NA | \(+\) | NA | NA | 1.749 | 0.154 |
| (Sum of weights) | 0.999 | 0.791 | 0.518 | 0.304 | 0.158 | 0.134 | NA | NA |
The average estimates for regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the 2 best meta-regression models. Note that mu shows the overall value at the intercept while beta represents the contrast (or slope) between two groups in the Unit column.
# getting averaged R2 and variance components not provided by the MuMIn package
average_sigma2 <- weighted.mean(x = c(model1$sigma2, model2$sigma2, model3$sigma2,
model4$sigma2), w = candidates_aic2$weight)
average_R2 <- weighted.mean(x = c(R2(model1)[1], R2(model2)[1], R2(model3)[1], R2(model4)[1]),
w = candidates_aic2$weight)
# creating a table
tibble(`Fixed effect` = c("Intercept (both-Microbe-Mutualist)", "Microbe-Plant",
"Microbe-Invert", "Microbe-Vert", "host_range", "both-horizontal", "both-vertical",
"Mutualist-Parasite"), Estimate = mr_averaged_aic2$coefmat.full[, 1], `Lower CI [0.025]` = mr_averaged_aic2$coefmat.full[,
1] - mr_averaged_aic2$coefmat.full[, 2] * qnorm(0.975), `Upper CI [0.975]` = mr_averaged_aic2$coefmat.full[,
1] + mr_averaged_aic2$coefmat.full[, 2] * qnorm(0.975), `V[authors]` = c(average_sigma2,
rep(NA, 7)), R2 = c(average_R2, rep(NA, 7))) %>% kable("html", digits = 3) %>%
kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| Intercept (both-Microbe-Mutualist) | 0.872 | 0.606 | 1.138 | 0.064 | 0.259 |
| Microbe-Plant | -0.415 | -0.685 | -0.144 | NA | NA |
| Microbe-Invert | -0.311 | -0.581 | -0.041 | NA | NA |
| Microbe-Vert | -0.278 | -0.533 | -0.023 | NA | NA |
| host_range | -0.077 | -0.215 | 0.060 | NA | NA |
| both-horizontal | 0.166 | -0.078 | 0.411 | NA | NA |
| both-vertical | -0.040 | -0.187 | 0.106 | NA | NA |
| Mutualist-Parasite | 0.014 | -0.078 | 0.106 | NA | NA |
Here, we conducted 3 kinds of publication bias analyses: 1) contour-enhanced funnel plots (Peters et al. 2008) of residuals (Egger et al. 1997; Nakagawa & Santos 2012), 2) a type of Egger regression (Egger et al. 1997; Moreno et al. 2009), and 3) a regression-based time-lag bias test (Nakagawa & Santos 2012).
A normal funnel plot assumes homogeneity (i.e., I2 = 0). Therefore, we controlled for important moderators (i.e., mode_of_transmission_broad, host_tax_broad, log(host_range_link_ratio), & symbiosis).
We do not observe normal skewness in our enhanced-counter funnel plot (Supplementary Figure 5). This funnel asymmetry seems different from one caused by publication bias(Peters et al. 2008); we do not expect a “hollow” in the region with high precision or i.e. inverse standard error (2.4-4) and relative high effect sizes (Zr = 0.5-1.0). The funnel asymmetry is mainly caused by the boundary created by the number of randomizations (see the “Sensitivity Analysis” section where we deal with this skewness).
#
res_funnel_plot <- rma.mv(yi = Zr, V = VZr, mods = ~mode_of_transmission_broad +
host_tax_broad + log(host_range_link_ratio) + symbiosis, random = ~1 | authors,
data = dat)
funnel(res_funnel_plot, yaxis = "seinv", level = c(90, 95, 99), xlim = c(-3, 4),
ylim = c(0.9, 5.2), xlab = "Residuals (correlation)", ylab = "Precision (1/SE)",
shade = c("white", "gray55", "gray75"), refline = 0, legend = TRUE)
Figure 5a: A residual funnel plot from the meta-regression model with mode_of_transmission_broad, host_tax_broad, & symbiosis; ‘residual value’ is on Zr and ‘inverse standard error’ is precision 1/sqrt(VZr).
Further, Egger regression analyses (see below) showed that sqrt(VZr) (sampling errors [SE] for effect sizes) accounts for much heterogeneity, so we added that to our model. The funnel asymmetry we see in Supplementary Figure 6 (if any) is much less severe than that in Supplementary Figure 5.
#
res_funnel_plot2 <- rma.mv(yi = Zr, V = VZr, mods = ~sqrt(VZr) + mode_of_transmission_broad +
host_tax_broad + log(host_range_link_ratio) + symbiosis, random = ~1 | authors,
data = dat)
funnel(res_funnel_plot2, yaxis = "seinv", level = c(90, 95, 99), xlim = c(-3, 4),
ylim = c(0.9, 5.2), xlab = "Residuals (correlation)", ylab = "Precision (1/SE)",
shade = c("white", "gray55", "gray75"), refline = 0, legend = TRUE)
Figure 5b: A residual funnel plot from the meta-regression model with sqrt(VZr), mode_of_transmission_broad, host_tax_broad, log(host_range_link_ratio), & symbiosis; ‘residual value’ is on Zr and ‘inverse standard error’ is precision 1/sqrt(VZr).
# building fig 3 using patchwork
par(mfrow = c(1, 2))
funnel(res_funnel_plot, yaxis = "seinv", level = c(90, 95, 99), xlim = c(-3, 4),
ylim = c(0.9, 5.2), xlab = "Residuals (correlation)", ylab = "Precision (1/SE)",
shade = c("white", "gray55", "gray75"), refline = 0, legend = TRUE)
mtext("a", side = 3, cex = 1.5, line = 1, adj = 0)
funnel(res_funnel_plot2, yaxis = "seinv", level = c(90, 95, 99), xlim = c(-3, 4),
ylim = c(0.9, 5.2), xlab = "Residuals (correlation)", ylab = "Precision (1/SE)",
shade = c("white", "gray55", "gray75"), refline = 0, legend = TRUE)
mtext("b", side = 3, cex = 1.5, line = 1, adj = 0)
We applied Egger’s regression to test whether the funnel asymmetries we observe in our funnel plots are statistically significant or not.
The test (or sqrt(VZr)) is significant. However, as mentioned above, this is due to the boundary created by the number of randomizations; this boundary can be seen in Figure S4.1 below.
#
egger_regression_uni <- rma.mv(yi = Zr, V = VZr, mods = ~sqrt(VZr), random = ~1 |
authors, data = dat)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with sqrt(VZr).
# getting marginal R2
r2_egger_regression_uni <- R2(egger_regression_uni)
# getting estimates: name does not work for slopes
res_egger_regression_uni <- get_est(egger_regression_uni, mod = "sqrt(VZr)")
# creating a table
tibble(`Fixed effect` = c("Intercept", "sqrt(VZr)"), Estimate = c(res_egger_regression_uni$estimate),
`Lower CI [0.025]` = c(res_egger_regression_uni$lowerCL), `Upper CI [0.975]` = c(res_egger_regression_uni$upperCL),
`V[authors]` = c(egger_regression_uni$sigma2, NA), R2 = c(r2_egger_regression_uni[1],
NA)) %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| Intercept | 0.178 | 0.054 | 0.301 | 0.061 | 0.397 |
| sqrt(VZr) | 1.389 | 0.944 | 1.833 | NA | NA |
pred_egger_regression_uni <- predict.rma(egger_regression_uni)
# plotting
fit_egger_regression_uni <- dat %>% mutate(ymin = pred_egger_regression_uni$ci.lb,
ymax = pred_egger_regression_uni$ci.ub, ymin2 = pred_egger_regression_uni$cr.lb,
ymax2 = pred_egger_regression_uni$cr.ub, pred = pred_egger_regression_uni$pred) %>%
ggplot(aes(x = sqrt(VZr), y = Zr, size = (1/VZr) + 3)) + geom_point(shape = 21,
fill = "grey90") + # geom_ribbon(aes(ymin = ymin, ymax = ymax), fill = '#0072B2') + # not quite sure
# why this does not work
geom_smooth(aes(y = ymin2), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25,
colour = "#0072B2") + geom_smooth(aes(y = ymax2), method = "loess", se = FALSE,
lty = "dotted", lwd = 0.25, colour = "#0072B2") + geom_smooth(aes(y = ymin),
method = "loess", se = FALSE, lty = "dotted", lwd = 0.25, colour = "#D55E00") +
geom_smooth(aes(y = ymax), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25,
colour = "#D55E00") + geom_smooth(aes(y = pred), method = "loess", se = FALSE,
lty = "dashed", lwd = 0.5, colour = "black") + ylim(-1, 2) + xlim(0.05, 0.45) +
# geom_abline(intercept = mr_host_range_link_ratio$beta[[1]], slope =
# mr_host_range_link_ratio$beta[[2]], alpha = 0.7, linetype = 'dashed', size =
# 0.5) +
labs(x = "sqrt(sampling variance)", y = expression(paste(italic(Zr), " (effect size)")),
size = expression(paste(italic(n), " (# of species pairs)"))) + guides(fill = "none",
colour = "none") + # themses
theme_bw() + theme(legend.position = c(0, 1), legend.justification = c(0, 1)) + theme(legend.direction = "horizontal") +
# theme(legend.background = element_rect(fill = 'white', colour = 'black')) +
theme(legend.background = element_blank()) + theme(axis.text.y = element_text(size = 10,
colour = "black", hjust = 0.5, angle = 90))
fit_egger_regression_uni
Figure S4.1: A bubble plot showing a predicted regression line for the continuous variable sqrt(VZr), indicating 95% confidence regions (orange dotted lines) and 95% prediction regions (blue dotted lines), with observed effect sizes based on various sample sizes.
We also conducted an Egger regression controlling other important moderators (i.e., mode_of_transmission_broad, host_tax_broad, log(host_range_link_ratio), & symbiosis). After controlling for these variables, sqrt(VZr) stays significant.
#
egger_regression_mul <- rma.mv(yi = Zr, V = VZr, mods = ~sqrt(VZr) + mode_of_transmission_broad +
host_tax_broad + log(host_range_link_ratio) + symbiosis, random = ~1 | authors,
data = dat)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with sqrt(VZr).
# getting marginal R2
r2_egger_regression_mul <- R2(egger_regression_mul)
# creating a table
tibble(`Fixed effect` = c("Intercept (both-Microbe-Mutualist)", "sqrt(VZr)", "both-horizontal",
"both-vertical", "Microbe-Plant", "Microbe-Invert", "Microbe-Vert", "host_range",
"Mutualist-Parasite"), Estimate = c(egger_regression_mul$b), `Lower CI [0.025]` = c(egger_regression_mul$ci.lb),
`Upper CI [0.975]` = c(egger_regression_mul$ci.ub), `V[authors]` = c(egger_regression_mul$sigma2,
rep(NA, 8)), R2 = c(r2_egger_regression_mul[1], rep(NA, 8))) %>% kable("html",
digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| Intercept (both-Microbe-Mutualist) | 0.420 | 0.123 | 0.717 | 0.044 | 0.583 |
| sqrt(VZr) | 1.344 | 0.908 | 1.779 | NA | NA |
| both-horizontal | -0.066 | -0.195 | 0.062 | NA | NA |
| both-vertical | 0.213 | -0.027 | 0.453 | NA | NA |
| Microbe-Plant | -0.287 | -0.546 | -0.028 | NA | NA |
| Microbe-Invert | -0.239 | -0.496 | 0.018 | NA | NA |
| Microbe-Vert | -0.146 | -0.391 | 0.100 | NA | NA |
| host_range | 0.052 | -0.082 | 0.186 | NA | NA |
| Mutualist-Parasite | -0.071 | -0.237 | 0.096 | NA | NA |
pred_egger_regression_mul <-predict.rma(egger_regression_mul)
# plotting
fit_egger_regression_mul <- dat %>%
filter(!is.na(mode_of_transmission_broad) & !is.na(host_tax_broad) & !is.na(symbiosis) & !is.na(host_range_link_ratio)) %>% # getting ride of NA values
mutate(ymin = pred_egger_regression_mul$ci.lb,
ymax = pred_egger_regression_mul$ci.ub,
ymin2 = pred_egger_regression_mul$cr.lb,
ymax2 = pred_egger_regression_mul$cr.ub,
pred = pred_egger_regression_mul$pred) %>%
ggplot(aes(x = sqrt(VZr), y = Zr, size = (1/VZr) + 3)) +
geom_point(shape = 21, fill = "grey90") +
#geom_ribbon(aes(ymin = ymin, ymax = ymax), fill = "#0072B2") + # not quite sure why this does not work
geom_smooth(aes(y = ymin2), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25, colour = "#0072B2") +
geom_smooth(aes(y = ymax2), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25, colour = "#0072B2") +
geom_smooth(aes(y = ymin), method = "loess", se = FALSE,lty = "dotted", lwd = 0.25, colour ="#D55E00") +
geom_smooth(aes(y = ymax), method = "loess", se = FALSE, lty ="dotted", lwd = 0.25, colour ="#D55E00") +
geom_smooth(aes(y = pred), method = "loess", se = FALSE, lty ="dashed", lwd = 0.5, colour ="black") +
ylim(-1, 2) + xlim(0.05, 0.45) +
#geom_abline(intercept = mr_host_range_link_ratio$beta[[1]], slope = mr_host_range_link_ratio$beta[[2]], alpha = 0.7, linetype = "dashed", size = 0.5) +
labs(x = "sqrt(sampling variance)", y = expression(paste(italic(Zr), " (effect size)")), size = expression(paste(italic(n), " (# of species pairs)"))) +
guides(fill = "none", colour = "none") +
# themses
theme_bw() +
theme(legend.position= c(0, 1), legend.justification = c(0, 1)) +
theme(legend.direction="horizontal") +
#theme(legend.background = element_rect(fill = "white", colour = "black")) +
theme(legend.background = element_blank()) +
theme(axis.text.y = element_text(size = 10, colour ="black", hjust = 0.5, angle = 90))
fit_egger_regression_mul
Figure S4.2: A bubble plot showing a predicted loess line for the continuous variable sqrt(VZr) (given the values of the other 3 variables in the model), with their 95% confidence regions (orange dotted lines) and 95% prediction regions (blue dotted lines) with observed effect sizes based on various sample sizes. Note that the lines are not linear as these are based on multivariate predictions of the data points.
We do not find any evidence of a time-lag effect (a decline in the magnitude of the effect over time) in either the univariate or multivariate models (Figure S4.X and Figure S4.X).
#
time_lag_effect_uni <- rma.mv(yi = Zr, V = VZr, mods = ~year, random = ~1 | authors,
data = dat)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with year.
# getting marginal R2
r2_time_lag_effect_uni <- R2(time_lag_effect_uni)
# getting estimates: name does not work for slopes
res_time_lag_effect_uni <- get_est(time_lag_effect_uni, mod = "year")
# creating a table
tibble(`Fixed effect` = c("Intercept", "Year"), Estimate = c(res_time_lag_effect_uni$estimate),
`Lower CI [0.025]` = c(res_time_lag_effect_uni$lowerCL), `Upper CI [0.975]` = c(res_time_lag_effect_uni$upperCL),
`V[authors]` = c(time_lag_effect_uni$sigma2, NA), R2 = c(r2_time_lag_effect_uni[1],
NA)) %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| Intercept | -3.732 | -23.253 | 15.789 | 0.082 | 0.002 |
| Year | 0.002 | -0.008 | 0.012 | NA | NA |
pred_time_lag_effect_uni <- predict.rma(time_lag_effect_uni)
# plotting
fit_time_lag_effect <- dat %>% mutate(ymin = pred_time_lag_effect_uni$ci.lb, ymax = pred_time_lag_effect_uni$ci.ub,
ymin2 = pred_time_lag_effect_uni$cr.lb, ymax2 = pred_time_lag_effect_uni$cr.ub,
pred = pred_time_lag_effect_uni$pred) %>% ggplot(aes(x = year, y = Zr, size = (1/VZr) +
3)) + geom_point(shape = 21, fill = "grey90") + # geom_ribbon(aes(ymin = ymin, ymax = ymax), fill = '#0072B2') + # not quite sure
# why this does not work
geom_smooth(aes(y = ymin2), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25,
colour = "#0072B2") + geom_smooth(aes(y = ymax2), method = "loess", se = FALSE,
lty = "dotted", lwd = 0.25, colour = "#0072B2") + geom_smooth(aes(y = ymin),
method = "loess", se = FALSE, lty = "dotted", lwd = 0.25, colour = "#D55E00") +
geom_smooth(aes(y = ymax), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25,
colour = "#D55E00") + geom_smooth(aes(y = pred), method = "loess", se = FALSE,
lty = "dashed", lwd = 0.5, colour = "black") + ylim(-1, 2) + xlim(1994, 2019) +
scale_x_continuous(breaks = c(1995, 2000, 2005, 2010, 2015, 2020)) + # geom_abline(intercept = mr_host_range_link_ratio$beta[[1]], slope =
# mr_host_range_link_ratio$beta[[2]], alpha = 0.7, linetype = 'dashed', size =
# 0.5) +
labs(x = "Year", y = expression(paste(italic(Zr), " (effect size)")), size = expression(paste(italic(n),
" (# of species pairs)"))) + guides(fill = "none", colour = "none") + # themses
theme_bw() + theme(legend.position = c(0, 1), legend.justification = c(0, 1)) + theme(legend.direction = "horizontal") +
# theme(legend.background = element_rect(fill = 'white', colour = 'black')) +
theme(legend.background = element_blank()) + theme(axis.text.y = element_text(size = 10,
colour = "black", hjust = 0.5, angle = 90))
fit_time_lag_effect
Figure S4.3: A bubble plot showing a predicted regression line for the contenious variable year, indicating 95% confidence regions (orange dotted lines) and 95% prediction regions (blue dotted lines), with observed effect sizes based on various sample sizes.
#
time_lag_effect_mul <- rma.mv(yi = Zr, V = VZr, mods = ~year + mode_of_transmission_broad +
host_tax_broad + log(host_range_link_ratio) + symbiosis, random = ~1 | authors,
data = dat)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the meta-regression with year.
# getting marginal R2
r2_time_lag_effect_mul <- R2(time_lag_effect_mul)
# creating a table
tibble(`Fixed effect` = c("Intercept (both-Microbe-Mutualist)", "Year", "both-horizontal",
"both-vertical", "Microbe-Plant", "Microbe-Invert", "Microbe-Vert", "host_range",
"Mutualist-Parasite"), Estimate = c(time_lag_effect_mul$b), `Lower CI [0.025]` = c(time_lag_effect_mul$ci.lb),
`Upper CI [0.975]` = c(time_lag_effect_mul$ci.ub), `V[authors]` = c(time_lag_effect_mul$sigma2,
rep(NA, 8)), R2 = c(r2_time_lag_effect_mul[1], rep(NA, 8))) %>% kable("html",
digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| Intercept (both-Microbe-Mutualist) | 2.773 | -16.237 | 21.783 | 0.066 | 0.259 |
| Year | -0.001 | -0.010 | 0.009 | NA | NA |
| both-horizontal | -0.073 | -0.217 | 0.071 | NA | NA |
| both-vertical | 0.154 | -0.109 | 0.417 | NA | NA |
| Microbe-Plant | -0.444 | -0.725 | -0.163 | NA | NA |
| Microbe-Invert | -0.322 | -0.601 | -0.044 | NA | NA |
| Microbe-Vert | -0.280 | -0.545 | -0.015 | NA | NA |
| host_range | 0.050 | -0.102 | 0.201 | NA | NA |
| Mutualist-Parasite | -0.071 | -0.253 | 0.110 | NA | NA |
pred_time_lag_effect_mul <- predict.rma(time_lag_effect_mul)
# plotting
fit_time_lag_effect_mul <- dat %>% filter(!is.na(mode_of_transmission_broad) & !is.na(host_tax_broad) &
!is.na(symbiosis) & !is.na(host_range_link_ratio)) %>% mutate(ymin = pred_time_lag_effect_mul$ci.lb,
ymax = pred_time_lag_effect_mul$ci.ub, ymin2 = pred_time_lag_effect_mul$cr.lb,
ymax2 = pred_time_lag_effect_mul$cr.ub, pred = pred_time_lag_effect_mul$pred) %>%
ggplot(aes(x = year, y = Zr, size = (1/VZr) + 3)) + geom_point(shape = 21, fill = "grey90") +
# geom_ribbon(aes(ymin = ymin, ymax = ymax), fill = '#0072B2') + # not quite sure
# why this does not work
geom_smooth(aes(y = ymin2), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25,
colour = "#0072B2") + geom_smooth(aes(y = ymax2), method = "loess", se = FALSE,
lty = "dotted", lwd = 0.25, colour = "#0072B2") + geom_smooth(aes(y = ymin),
method = "loess", se = FALSE, lty = "dotted", lwd = 0.25, colour = "#D55E00") +
geom_smooth(aes(y = ymax), method = "loess", se = FALSE, lty = "dotted", lwd = 0.25,
colour = "#D55E00") + geom_smooth(aes(y = pred), method = "loess", se = FALSE,
lty = "dashed", lwd = 0.5, colour = "black") + ylim(-1, 2) + xlim(1994, 2019) +
scale_x_continuous(breaks = c(1995, 2000, 2005, 2010, 2015, 2020)) + # geom_abline(intercept = mr_host_range_link_ratio$beta[[1]], slope =
# mr_host_range_link_ratio$beta[[2]], alpha = 0.7, linetype = 'dashed', size =
# 0.5) +
labs(x = "Year", y = expression(paste(italic(Zr), " (effect size)")), size = expression(paste(italic(n),
" (# of species pairs)"))) + guides(fill = "none", colour = "none") + # themses
theme_bw() + theme(legend.position = c(0, 1), legend.justification = c(0, 1)) + theme(legend.direction = "horizontal") +
# theme(legend.background = element_rect(fill = 'white', colour = 'black')) +
theme(legend.background = element_blank()) + theme(axis.text.y = element_text(size = 10,
colour = "black", hjust = 0.5, angle = 90))
fit_time_lag_effect_mul
Figure S4.4: A bubble plot showing a predicted loess line for the continuous variable year (given the values of the other 3 variables in the model), indicating 95% confidence regions (orange dotted lines) and 95% prediction regions (blue dotted lines) with observed effect sizes based on various sample sizes. Note that the lines are not linear as these are based on multivariate predictions of the data points.
The funnel plots above identified the issue of upper bounds for the effect size given a sample size (an upper limit of a p value given the number of randomizations). This boundary would influence our estimates of mean effect sizes and contrasts (i.e., comparing two groups), often making our overall conclusions too conservative. To demonstrate this, we conducted two analyses to show: 1) the number of randomizations (log(no_randomizations)) do not differ between categories in the 3 important categorical moderators (mode_of_transmission_broad, host_tax_broad, & symbiosis), and, 2) categories with high effect sizes would include “bounded” effect sizes (i.e., from p = 0.01, 0.001, or 0.0001; limit_reached) in the 3 moderators.
Below, we showed that none of the categories have significantly different numbers of randomizations in all mode_of_transmission_broad, host_tax_broad, & symbiosis.
# 233 --- Yes = 74 (0.3175966%); No = 159
# symbiosis multiple intercepts
sa_random_symbiosis1 <- lmer(log(no_randomizations) ~ symbiosis - 1 + (1 | authors),
data = dat)
# contrast
sa_random_symbiosis2 <- lmer(log(no_randomizations) ~ symbiosis + (1 | authors),
data = dat)
Regression coefficients (Estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2), from the regression with symbiosis on log(no_randomizations).
# getting marginal R2
r2_sa_random_symbiosis <- r2_nakagawa(sa_random_symbiosis1)
# getting estimates
res_sa_random_symbiosis <- tibble(estiamte = c(fixef(sa_random_symbiosis1), fixef(sa_random_symbiosis2)[2]))
ci_sa_random_symbiosis1 <- confint(sa_random_symbiosis1)
ci_sa_random_symbiosis2 <- confint(sa_random_symbiosis2)
res_sa_random_symbiosis %<>% mutate(lowerCL = c(ci_sa_random_symbiosis1[3:4, 1],
ci_sa_random_symbiosis2[4, 1]))
res_sa_random_symbiosis %<>% mutate(upperCL = c(ci_sa_random_symbiosis1[3:4, 2],
ci_sa_random_symbiosis2[4, 2]))
# creating a table
tibble(`Fixed effect` = c(as.character(res_symbiosis1$name), cont_gen(res_symbiosis1$name)),
Estimate = res_sa_random_symbiosis$estiamte, `Lower CI [0.025]` = res_sa_random_symbiosis$lowerCL,
`Upper CI [0.975]` = res_sa_random_symbiosis$upperCL, `V[authors]` = c(attr(VarCorr(sa_random_symbiosis1)$author,
"stddev")^2, rep(NA, 2)), `V[residuals]` = c(attr(VarCorr(sa_random_symbiosis1),
"sc")^2, rep(NA, 2)), R2 = c(r2_sa_random_symbiosis$R2_marginal, rep(NA,
2))) %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | V[residuals] | R2 |
|---|---|---|---|---|---|---|
| Mutualist | 7.696 | 7.286 | 8.106 | 2.931 | 0.352 | 0.004 |
| Parasite | 7.619 | 7.314 | 7.924 | NA | NA | NA |
| Mutualist-Parasite | -0.077 | -0.544 | 0.389 | NA | NA | NA |
# host_tax_broad mutiple intercepts
sa_random_host_tax_broad1 <- lmer(log(no_randomizations) ~ host_tax_broad - 1 + (1 |
authors), data = dat)
# contrast 1
sa_random_host_tax_broad2 <- lmer(log(no_randomizations) ~ host_tax_broad + (1 |
authors), data = dat)
# contrast 2
sa_random_host_tax_broad3 <- lmer(log(no_randomizations) ~ relevel(host_tax_broad,
ref = "Plant") + (1 | authors), data = dat)
# contrast 3
sa_random_host_tax_broad4 <- lmer(log(no_randomizations) ~ relevel(host_tax_broad,
ref = "Invert") + (1 | authors), data = dat)
Regression coefficients (estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the regression with host_tax_broad on log(no_randomizations).
# getting marginal R2
r2_sa_random_host_tax_broad <- r2_nakagawa(sa_random_host_tax_broad1)
# getting estimates
res_sa_random_host_tax_broad <- tibble(estiamte = c(fixef(sa_random_host_tax_broad1),
fixef(sa_random_host_tax_broad2)[2:4],
fixef(sa_random_host_tax_broad3)[3:4],
fixef(sa_random_host_tax_broad4)[4]))
ci_sa_random_host_tax_broad1<-confint(sa_random_host_tax_broad1)
ci_sa_random_host_tax_broad2<-confint(sa_random_host_tax_broad2)
ci_sa_random_host_tax_broad3<-confint(sa_random_host_tax_broad3)
ci_sa_random_host_tax_broad4<-confint(sa_random_host_tax_broad4)
res_sa_random_host_tax_broad %<>% mutate(lowerCL = c(ci_sa_random_host_tax_broad1[3:6,1],
ci_sa_random_host_tax_broad2[4:6,1],
ci_sa_random_host_tax_broad3[5:6,1],
ci_sa_random_host_tax_broad4[6,1]))
res_sa_random_host_tax_broad %<>% mutate(upperCL = c(ci_sa_random_host_tax_broad1[3:6,2],
ci_sa_random_host_tax_broad2[4:6,2],
ci_sa_random_host_tax_broad3[5:6,2],
ci_sa_random_host_tax_broad4[6,2]))
# creating a table
tibble(
`Fixed effect` = c(as.character(res_symbiont_tax_broad1$name), cont_gen(res_symbiont_tax_broad1$name)), # done
Estimate = res_sa_random_host_tax_broad$estiamte,
`Lower CI [0.025]` = res_sa_random_host_tax_broad$lowerCL,
`Upper CI [0.975]` = res_sa_random_host_tax_broad$upperCL,
`V[authors]` = c(attr(VarCorr(sa_random_host_tax_broad1)$author,"stddev")^2, rep(NA, 9)),
`V[residuals]` =c(attr(VarCorr(sa_random_host_tax_broad1),"sc")^2, rep(NA, 9)),
`R2` = c(r2_sa_random_host_tax_broad$R2_marginal, rep(NA, 9))) %>% kable("html", digits = 3) %>%
kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | V[residuals] | R2 |
|---|---|---|---|---|---|---|
| Microbe | 8.223 | 7.279 | 9.168 | 2.937 | 0.35 | 0.071 |
| Plant | 7.466 | 6.856 | 8.076 | NA | NA | NA |
| Invert | 7.650 | 7.167 | 8.134 | NA | NA | NA |
| Vert | 7.634 | 7.265 | 8.003 | NA | NA | NA |
| Microbe-Plant | -0.757 | -1.882 | 0.367 | NA | NA | NA |
| Microbe-Invert | -0.573 | -1.634 | 0.488 | NA | NA | NA |
| Microbe-Vert | -0.589 | -1.603 | 0.425 | NA | NA | NA |
| Plant-Invert | 0.184 | -0.594 | 0.962 | NA | NA | NA |
| Plant-Vert | 0.168 | -0.545 | 0.881 | NA | NA | NA |
| Invert-Vert | -0.016 | -0.599 | 0.566 | NA | NA | NA |
# mode_of_transmission_broad
sa_random_mode_of_transmission_broad1 <- lmer(log(no_randomizations) ~ mode_of_transmission_broad -
1 + (1 | authors), data = dat)
# contrast 1
sa_random_mode_of_transmission_broad2 <- lmer(log(no_randomizations) ~ mode_of_transmission_broad +
(1 | authors), data = dat)
# contrast 2
sa_random_mode_of_transmission_broad3 <- lmer(log(no_randomizations) ~ relevel(mode_of_transmission_broad,
ref = "vertical") + (1 | authors), data = dat)
Regression coefficients (estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the regression with mode_of_transmission_broad on log(no_randomizations).
# getting marginal R2
r2_sa_random_mode_of_transmission_broad <- r2_nakagawa(sa_random_mode_of_transmission_broad1)
# getting estimates
res_sa_random_mode_of_transmission_broad <- tibble(estiamte = c(fixef(sa_random_mode_of_transmission_broad1),
fixef(sa_random_mode_of_transmission_broad2)[2:3], fixef(sa_random_mode_of_transmission_broad3)[3]))
ci_sa_random_mode_of_transmission_broad1 <- confint(sa_random_mode_of_transmission_broad1)
ci_sa_random_mode_of_transmission_broad2 <- confint(sa_random_mode_of_transmission_broad2)
ci_sa_random_mode_of_transmission_broad3 <- confint(sa_random_mode_of_transmission_broad3)
res_sa_random_mode_of_transmission_broad %<>% mutate(lowerCL = c(ci_sa_random_mode_of_transmission_broad1[3:5,
1], ci_sa_random_mode_of_transmission_broad2[4:5, 1], ci_sa_random_mode_of_transmission_broad3[5,
1]))
res_sa_random_mode_of_transmission_broad %<>% mutate(upperCL = c(ci_sa_random_mode_of_transmission_broad1[3:5,
2], ci_sa_random_mode_of_transmission_broad2[4:5, 2], ci_sa_random_mode_of_transmission_broad3[5,
2]))
# creating a table
tibble(`Fixed effect` = c(as.character(res_mode_of_transmission_broad1$name), cont_gen(res_mode_of_transmission_broad1$name)),
Estimate = res_sa_random_mode_of_transmission_broad$estiamte, `Lower CI [0.025]` = res_sa_random_mode_of_transmission_broad$lowerCL,
`Upper CI [0.975]` = res_sa_random_mode_of_transmission_broad$upperCL, `V[authors]` = c(attr(VarCorr(sa_random_mode_of_transmission_broad1)$author,
"stddev")^2, rep(NA, 5)), `V[residuals]` = c(attr(VarCorr(sa_random_mode_of_transmission_broad1),
"sc")^2, rep(NA, 5)), R2 = c(r2_sa_random_mode_of_transmission_broad$R2_marginal,
rep(NA, 5))) %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | V[residuals] | R2 |
|---|---|---|---|---|---|---|
| both | 7.535 | 6.978 | 8.091 | 2.913 | 0.363 | 0.12 |
| horizontal | 7.515 | 7.164 | 7.866 | NA | NA | NA |
| vertical | 8.089 | 7.509 | 8.668 | NA | NA | NA |
| both-horizontal | -0.020 | -0.678 | 0.638 | NA | NA | NA |
| both-vertical | 0.554 | -0.250 | 1.357 | NA | NA | NA |
| horizontal-vertical | -0.574 | -1.251 | 0.104 | NA | NA | NA |
Below, we show that categories with higher effect sizes were more likely to have “bounded” effect sizes (limit_reached) in all mode_of_transmission_broad, host_tax_broad, & symbiosis. This indicates that the higher the estimate of mean effect size is, the more underestimated the mean effect size is. This is true for differences between two categories; the larger the difference between the two, the more underestimated the difference is.
Regression coefficients (estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the regression with symbiosis on limit_reached (compare the effect estimates in Table S3.1 with the corresponding values in this table below; their ranking in estimation should usually match).
# symbiosis
sa_limit_symbiosis1 <- glmer(limit_reached ~ symbiosis - 1 + (1 | authors), family = "binomial",
data = dat)
# getting marginal R2
r2_sa_limit_symbiosis <- r2_nakagawa(sa_limit_symbiosis1)
# getting estimates
res_sa_limit_symbiosis <- tibble(estiamte = fixef(sa_limit_symbiosis1))
res_sa_limit_symbiosis %<>% mutate(lowerCL = (tidy(sa_limit_symbiosis1)$estimate[-3] -
tidy(sa_limit_symbiosis1)$std.error[-3] * qnorm(0.975)))
res_sa_limit_symbiosis %<>% mutate(upperCL = (tidy(sa_limit_symbiosis1)$estimate[-3] +
tidy(sa_limit_symbiosis1)$std.error[-3] * qnorm(0.975)))
# creating a table
tibble(`Fixed effect` = as.character(res_symbiosis1$name), Estimate = res_sa_limit_symbiosis$estiamte,
`Lower CI [0.025]` = res_sa_limit_symbiosis$lowerCL, `Upper CI [0.975]` = res_sa_limit_symbiosis$upperCL,
`V[authors]` = c(attr(VarCorr(sa_limit_symbiosis1)$author, "stddev")^2, rep(NA,
1)), R2 = c(r2_sa_limit_symbiosis$R2_marginal, rep(NA, 1))) %>% kable("html",
digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| Mutualist | -0.401 | -1.074 | 0.272 | 1.657 | 0.051 |
| Parasite | -1.309 | -2.015 | -0.603 | NA | NA |
# t_symbiosis
In both tables, mutualist has higher estimates (i.e, mutualist rearched the ceiling more often)
t_symbiosis %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| Mutualist | Zr (mu) | 0.631 | 0.534 | 0.727 | 0.078 | 0.054 |
| Parasite | Zr (mu) | 0.487 | 0.419 | 0.554 | NA | NA |
| Mutualist | r (mu) | 0.559 | 0.488 | 0.621 | NA | NA |
| Parasite | r (mu) | 0.451 | 0.396 | 0.503 | NA | NA |
| Mutualist-Parasite | Zr (beta) | -0.144 | -0.260 | -0.028 | NA | NA |
Regression coefficients (estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the regression with host_tax_broad on limit_reached (compare the effect estimates in Table S3.2 with the corresponding values in this table below; their ranking in estimation should usually match).
# host_tax_broad
sa_limit_host_tax_broad1 <- glmer(limit_reached ~ host_tax_broad - 1 + (1 | authors),
family = "binomial", data = dat)
# getting marginal R2
r2_sa_limit_host_tax_broad <- r2_nakagawa(sa_limit_host_tax_broad1)
# getting estimates
res_sa_limit_host_tax_broad <- tibble(estiamte = fixef(sa_limit_host_tax_broad1))
res_sa_limit_host_tax_broad %<>% mutate(lowerCL = (tidy(sa_limit_host_tax_broad1)$estimate[-5] -
tidy(sa_limit_host_tax_broad1)$std.error[-5] * qnorm(0.975)))
res_sa_limit_host_tax_broad %<>% mutate(upperCL = (tidy(sa_limit_host_tax_broad1)$estimate[-5] +
tidy(sa_limit_host_tax_broad1)$std.error[-5] * qnorm(0.975)))
# creating a table
tibble(`Fixed effect` = as.character(res_symbiont_tax_broad1$name), Estimate = res_sa_limit_host_tax_broad$estiamte,
`Lower CI [0.025]` = res_sa_limit_host_tax_broad$lowerCL, `Upper CI [0.975]` = res_sa_limit_host_tax_broad$upperCL,
`V[authors]` = c(attr(VarCorr(sa_limit_host_tax_broad1)$author, "stddev")^2,
rep(NA, 3)), R2 = c(r2_sa_limit_host_tax_broad$R2_marginal, rep(NA, 3))) %>%
kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| Microbe | -1.211 | -2.737 | 0.316 | 1.358 | 0.035 |
| Plant | -1.574 | -2.571 | -0.578 | NA | NA |
| Invert | -0.579 | -1.308 | 0.150 | NA | NA |
| Vert | -0.937 | -1.588 | -0.286 | NA | NA |
# t_host_tax
t_host_tax %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| Microbe | Zr (mu) | 0.883 | 0.639 | 1.126 | 0.077 | 0.162 |
| Plant | Zr (mu) | 0.377 | 0.256 | 0.498 | NA | NA |
| Invert | Zr (mu) | 0.629 | 0.522 | 0.736 | NA | NA |
| Vert | Zr (mu) | 0.515 | 0.435 | 0.595 | NA | NA |
| Microbe | r (mu) | 0.708 | 0.564 | 0.810 | NA | NA |
| Plant | r (mu) | 0.360 | 0.250 | 0.460 | NA | NA |
| Invert | r (mu) | 0.557 | 0.479 | 0.627 | NA | NA |
| Vert | r (mu) | 0.474 | 0.410 | 0.534 | NA | NA |
| Microbe-Plant | Zr (beta) | -0.506 | -0.778 | -0.234 | NA | NA |
| Microbe-Invert | Zr (beta) | -0.254 | -0.519 | 0.012 | NA | NA |
| Microbe-Vert | Zr (beta) | -0.367 | -0.624 | -0.111 | NA | NA |
| Plant-Invert | Zr (beta) | 0.252 | 0.091 | 0.414 | NA | NA |
| Plant-Vert | Zr (beta) | 0.139 | -0.006 | 0.284 | NA | NA |
| Invert-Vert | Zr (beta) | -0.114 | -0.244 | 0.017 | NA | NA |
Regression coefficients (estimate), 95% confidence intervals (CIs), variance components (V) and variance explained, R2[marginal] (R2) from the regression with mode_of_transmission_broad on limit_reached (compare the effect estimates in Table S3.7 with the corresponding values in this table below; their ranking in estimation should usually match).
# mode_of_transmission_broad
sa_limit_mode_of_transmission_broad1 <- glmer(limit_reached ~ mode_of_transmission_broad -
1 + (1 | authors), family = "binomial", data = dat)
# getting marginal R2
r2_sa_limit_mode_of_transmission_broad <- r2_nakagawa(sa_limit_mode_of_transmission_broad1)
# getting estimates
res_sa_limit_mode_of_transmission_broad <- tibble(estiamte = fixef(sa_limit_mode_of_transmission_broad1))
res_sa_limit_mode_of_transmission_broad %<>% mutate(lowerCL = (tidy(sa_limit_mode_of_transmission_broad1)$estimate[-4] -
tidy(sa_limit_mode_of_transmission_broad1)$std.error[-4] * qnorm(0.975)))
res_sa_limit_mode_of_transmission_broad %<>% mutate(upperCL = (tidy(sa_limit_mode_of_transmission_broad1)$estimate[-4] +
tidy(sa_limit_mode_of_transmission_broad1)$std.error[-4] * qnorm(0.975)))
# creating a table
tibble(`Fixed effect` = as.character(res_mode_of_transmission_broad1$name), Estimate = res_sa_limit_mode_of_transmission_broad$estiamte,
`Lower CI [0.025]` = res_sa_limit_mode_of_transmission_broad$lowerCL, `Upper CI [0.975]` = res_sa_limit_mode_of_transmission_broad$upperCL,
`V[authors]` = c(attr(VarCorr(sa_limit_mode_of_transmission_broad1)$author, "stddev")^2,
rep(NA, 2)), R2 = c(r2_sa_limit_mode_of_transmission_broad$R2_marginal, rep(NA,
2))) %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|
| both | -0.928 | -1.809 | -0.046 | 1.427 | 0.078 |
| horizontal | -1.326 | -2.031 | -0.621 | NA | NA |
| vertical | 0.060 | -0.750 | 0.869 | NA | NA |
# t_transmission
t_transmission %>% kable("html", digits = 3) %>% kable_styling("striped", position = "left")
| Fixed effect | Unit | Estimate | Lower CI [0.025] | Upper CI [0.975] | V[authors] | R2 |
|---|---|---|---|---|---|---|
| both | Zr (mu) | 0.578 | 0.463 | 0.693 | 0.061 | 0.187 |
| horizontal | Zr (mu) | 0.446 | 0.378 | 0.515 | NA | NA |
| vertical | Zr (mu) | 0.751 | 0.634 | 0.868 | NA | NA |
| both | r (mu) | 0.521 | 0.432 | 0.600 | NA | NA |
| horizontal | r (mu) | 0.419 | 0.361 | 0.474 | NA | NA |
| vertical | r (mu) | 0.636 | 0.561 | 0.701 | NA | NA |
| both-horizontal | Zr (beta) | -0.131 | -0.265 | 0.002 | NA | NA |
| both-vertical | Zr (beta) | 0.174 | 0.009 | 0.338 | NA | NA |
| horizontal-vertical | Zr (beta) | -0.305 | -0.441 | -0.169 | NA | NA |
Many coding materials have been borrowed from these papers (Cally et al. 2019; O’Dea et al. 2019). We thank Losia Lagisz for preparing small icons and cartoons used in the figures.
# pander for making it look nicer
sessionInfo() %>% pander()
R version 4.0.3 (2020-10-10)
Platform: x86_64-apple-darwin17.0 (64-bit)
locale: en_AU.UTF-8||en_AU.UTF-8||en_AU.UTF-8||C||en_AU.UTF-8||en_AU.UTF-8
attached base packages: grid, stats, graphics, grDevices, utils, datasets, methods and base
other attached packages: here(v.1.0.1), patchwork(v.1.1.1), png(v.0.1-7), performance(v.0.6.1), broom.mixed(v.0.2.6), lme4(v.1.1-26), MuMIn(v.1.43.17), plotly(v.4.9.3), ggbeeswarm(v.0.6.0), MCMCglmm(v.2.30), ape(v.5.4-1), coda(v.0.19-4), metafor(v.2.5-74), Matrix(v.1.3-2), pander(v.0.6.3), magrittr(v.2.0.1), gridExtra(v.2.3), kableExtra(v.1.3.1), forcats(v.0.5.1), stringr(v.1.4.0), dplyr(v.1.0.3), purrr(v.0.3.4), readr(v.1.4.0), tidyr(v.1.1.2), tibble(v.3.0.5), ggplot2(v.3.3.3) and tidyverse(v.1.3.0)
loaded via a namespace (and not attached): nlme(v.3.1-151), fs(v.1.5.0), lubridate(v.1.7.9.2), insight(v.0.12.0), webshot(v.0.5.2), httr(v.1.4.2), rprojroot(v.2.0.2), tensorA(v.0.36.2), tools(v.4.0.3), TMB(v.1.7.18), backports(v.1.2.1), R6(v.2.5.0), vipor(v.0.4.5), mgcv(v.1.8-33), DBI(v.1.1.1), lazyeval(v.0.2.2), colorspace(v.2.0-0), withr(v.2.4.1), tidyselect(v.1.1.0), compiler(v.4.0.3), cli(v.2.2.0), rvest(v.0.3.6), formatR(v.1.7), pacman(v.0.5.1), xml2(v.1.3.2), labeling(v.0.4.2), bayestestR(v.0.8.2), scales(v.1.1.1), digest(v.0.6.27), minqa(v.1.2.4), rmarkdown(v.2.6), pkgconfig(v.2.0.3), htmltools(v.0.5.1.1), highr(v.0.8), dbplyr(v.2.0.0), htmlwidgets(v.1.5.3), rlang(v.0.4.10), readxl(v.1.3.1), rstudioapi(v.0.13), farver(v.2.0.3), generics(v.0.1.0), jsonlite(v.1.7.2), Rcpp(v.1.0.6), munsell(v.0.5.0), fansi(v.0.4.2), lifecycle(v.0.2.0), stringi(v.1.5.3), yaml(v.2.2.1), mathjaxr(v.1.0-1), MASS(v.7.3-53), plyr(v.1.8.6), parallel(v.4.0.3), crayon(v.1.3.4), lattice(v.0.20-41), haven(v.2.3.1), splines(v.4.0.3), hms(v.1.0.0), knitr(v.1.31), pillar(v.1.4.7), boot(v.1.3-26), cubature(v.2.0.4.1), corpcor(v.1.6.9), codetools(v.0.2-18), reshape2(v.1.4.4), stats4(v.4.0.3), reprex(v.1.0.0), glue(v.1.4.2), evaluate(v.0.14), data.table(v.1.13.6), modelr(v.0.1.8), vctrs(v.0.3.6), nloptr(v.1.2.2.2), cellranger(v.1.1.0), gtable(v.0.3.0), assertthat(v.0.2.1), xfun(v.0.20), broom(v.0.7.3), viridisLite(v.0.3.0), beeswarm(v.0.2.3), statmod(v.1.4.35) and ellipsis(v.0.3.1)